Xyce  6.1
N_DEV_MemristorTEAM.C
Go to the documentation of this file.
1 // ----------------------------------------------------------------------------
2 // Copyright Notice
3 //
4 // Copyright 2002 Sandia Corporation. Under the terms
5 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S.
6 // Government retains certain rights in this software.
7 //
8 // Xyce(TM) Parallel Electrical Simulator
9 // Copyright (C) 2002-2015 Sandia Corporation
10 //
11 // This program is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // This program is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU General Public License
22 // along with this program. If not, see <http://www.gnu.org/licenses/>.
23 // ----------------------------------------------------------------------------
24 
25 //----------------------------------------------------------------------------
26 // Filename : $RCSfile: N_DEV_MemristorTEAM.C,v $
27 //
28 // Purpose : Implementation of the TEAM memristor model. See.
29 // TEAM: Threshold Adaptive Memristor Model
30 // Shahar Kvatinsky, Eby G. Friedman, Uri Weiser,
31 // IEEE Transactions on Circuits and Systems-I Vol 60, No. 1, Jan. 2013
32 // model details.
33 //
34 // Special Notes :
35 //
36 // Creator : Richard Schiek, Electrical Models & Simulations
37 //
38 // Creation Date : 10/23/14
39 //
40 // Revision Information:
41 // ---------------------
42 //
43 // Revision Number: $Revision: 1.14 $
44 //
45 // Revision Date : $Date: 2015/08/06 19:59:02 $
46 //
47 // Current Owner : $Author: rlschie $
48 //----------------------------------------------------------------------------
49 #include <Xyce_config.h>
50 
51 #include <N_DEV_MemristorTEAM.h>
52 
53 #include <N_DEV_DeviceOptions.h>
54 #include <N_DEV_ExternData.h>
55 #include <N_DEV_SolverState.h>
56 #include <N_DEV_Message.h>
57 #include <N_LAS_Matrix.h>
58 #include <N_UTL_FeatureTest.h>
59 #include <Sacado.hpp>
60 
61 namespace Xyce {
62 namespace Device {
63 namespace MemristorTEAM {
64 
65 //
66 // Window functions used to keep the internal variable
67 // "w" when it characterized the width of the conductive layer
68 // "x" when it characterizes the thickness of the conductive layer
69 //
70 //
71 template <typename ScalarT>
72 void JogelkarWindowFunction( const ScalarT & w, const ScalarT & D, const ScalarT & p, ScalarT & fval )
73 {
74  fval = 1.0 - pow( (2.0*(w/D) - 1.0), (2*p) );
75 }
76 
77 
78 template <typename ScalarT>
79 void BiolekWindowFunction( const ScalarT & w, const ScalarT & D, const ScalarT & p, const ScalarT & i, ScalarT & fval )
80 {
81  ScalarT stp_of_i = 1.0;
82  if( i < 0.0 )
83  stp_of_i = 0;
84  fval = 1.0 - pow( ((w/D) - stp_of_i), (2*p) );
85 }
86 
87 template <typename ScalarT>
88 void ProdromakisWindowFunction( const ScalarT & w, const ScalarT & D, const ScalarT & p, const ScalarT & j, ScalarT & fval )
89 {
90  fval = j * (1.0 - pow( (pow( (w - 0.5), 2) -0.75), p ) );
91 }
92 
93 
94 //
95 // Window function F for TEAM window function
96 // takes x, Aoff and wc as parameters.
97 // returns fval as result
98 template <typename ScalarT>
99 void TEAMWindowFunctionF(const ScalarT & x, const ScalarT & i, const ScalarT & aOff, ScalarT & aOn, const ScalarT & wc, ScalarT & fval )
100 {
101  if( i >= 0.0 )
102  {
103  fval = exp( -exp( (x-aOff)/wc ) );
104  }
105  else
106  {
107  fval = exp( -exp( -(x-aOn)/wc ) );
108  }
109 }
110 
111 
112 // Internal state variable equations needed by the TEAM model
113 //
114 // dx(t)/dt + F(X) = 0
115 //
116 // F(x) : defined as.
117 //
118 // F(x) = -k_off ( (i(t)/i_off) - 1)^a_off F_off(x) if 0 < i_off < i
119 //
120 // F(x) = 0 if i_on < i < i_off
121 //
122 // F(x) = -k_on ( (i(t)/i_on) - 1 )^a_on F_on(x) if i < i_on < 0
123 //
124 // Note F(x) has a "-" sign in the above relative to the paper it's
125 // based on because I've moved it to the left hand side of the equation.
126 
127 //
128 // State variable x. F vector term
129 //
130 // F(x) : defined as.
131 //
132 // F(x) = -k_off ( (i(t)/i_off) - 1)^a_off F_off(x) if 0 < i_off < i
133 //
134 // F(x) = 0 if i_on < i < i_off
135 //
136 // F(x) = -k_on ( (i(t)/i_on) - 1 )^a_on F_on(x) if i < i_on < 0
137 //
138 // Note F(x) has a "-" sign in the above relative to the paper it's
139 // based on because I've moved it to the left hand side of the equation.
140 
141 template <typename ScalarT>
142 void xVarFterm( const ScalarT & Vpos, const ScalarT & Vneg, const ScalarT & x,
143  const ScalarT & G, const ScalarT & iOff, const ScalarT & iOn,
144  const ScalarT & kOff, const ScalarT & kOn,
145  const ScalarT & alphaOff, const ScalarT & alphaOn, ScalarT & fval )
146 {
147 
148  // equation 14a, 14b and 14c
149  ScalarT i = G * (Vpos - Vneg);
150 
151  if( i>=iOff )
152  {
153  fval=-kOff*pow( ((i/iOff) - 1.0), alphaOff );
154  }
155  else if( i<=iOn )
156  {
157  fval=-kOn*pow( ((i/iOn) - 1.0), alphaOn);
158  }
159  if ( (i>iOn) && (i<iOff) )
160  {
161  fval=0.0;
162  }
163 }
164 
165 
166 // linear current voltage model, Reff
167 template <typename ScalarT>
168 void ReffLin( const ScalarT & X, const ScalarT & Ron, const ScalarT & Roff, const ScalarT & Xon, const ScalarT & Xoff, ScalarT & fval )
169 {
170  // equation 15
171  fval = Ron + (X - Xon) * (Roff-Ron)/(Xoff - Xon);
172 }
173 
174 // nonlinear current voltage model, Reff
175 template <typename ScalarT>
176 void ReffNonLin( const ScalarT & X, const ScalarT & Ron, const ScalarT & Roff, const ScalarT & Xon, const ScalarT & Xoff, ScalarT & fval )
177 {
178  // equation 17
179  ScalarT lambda = log(Roff/Ron);
180  // equation 16
181  fval = Ron * exp( (X - Xon) * lambda /(Xoff - Xon) );
182 }
183 
184 
185 // Common Jacobian Stamp for all MemristorTEAM devices.
186 // Because all memristors have identical Jacobian stamps, this data is
187 // declared static and is shared by all memristor instances.
188 //
189 std::vector<std::vector<int> > Instance::jacStamp;
190 
191 //-----------------------------------------------------------------------------
192 // Function : Xyce::Device::MemristorTEAM::Instance::initializeJacobianStamp
193 // Purpose :
194 // Special Notes :
195 // Scope : private
196 // Creator : Richard Schiek, Electrical Models & Simulations
197 // Creation Date : 2/11/2014
198 //-----------------------------------------------------------------------------
199 //
200 // @brief Common Jacobian stamp initializer for all MemristorTEAM devices.
201 //
202 // The Jacobian stamp is a sparse-matrix representation of the pattern
203 // of non-zero elements that a device will put into the Jacobian matrix.
204 //
205 // The Jacobian stamp is used by the Topology package to determine indices
206 // into the full Jacobian matrix for elements that correspond to this
207 // device.
208 //
209 // There is one row of the Jacobian stamp for each equation associated with
210 // a device. The number of elements in a row are the number of non-zero
211 // elements in that row of the device's contribution to the Jacobian.
212 // The values of the elements are numbers local to the device that
213 // represent the column in which the non-zero element appears.
214 //
215 // For this memristor, there are two external nodes (the positive and negative
216 // terminals of the device). The positive node is referred to as the 0th
217 // node of the device, and the negative node the 1st node.
218 //
220 {
221  if (jacStamp.empty())
222  {
223  jacStamp.resize(3);
224  jacStamp[0].resize(3);
225  jacStamp[0][0] = 0;
226  jacStamp[0][1] = 1;
227  jacStamp[0][2] = 2;
228  jacStamp[1].resize(3);
229  jacStamp[1][0] = 0;
230  jacStamp[1][1] = 1;
231  jacStamp[1][2] = 2;
232  jacStamp[2].resize(3);
233  jacStamp[2][0] = 0;
234  jacStamp[2][1] = 1;
235  jacStamp[2][2] = 2;
236  }
237 }
238 
239 //-----------------------------------------------------------------------------
240 // Function : Xyce::Device::MemristorTEAM::Traits::loadInstanceParameters
241 // Purpose :
242 // Special Notes :
243 // Scope : private
244 // Creator : Richard Schiek, Electrical Models & Simulations
245 // Creation Date : 10/23/2014
246 //-----------------------------------------------------------------------------
247 //
248 // Loads the parameter definition into the instance parameter map.
249 //
250 // @param p instance parameter map
251 //
252 // Each parameter supported by the memristor device instance is
253 // defined via the addPar call. The minimum information required is
254 // the name of the parameter, the default value, and a member pointer
255 // to a variable in the instance class to which the value will be
256 // assigned.
257 //
258 // Additional information such as documentation for the parameter, units
259 // (used only in output of tables for documentation), whether a special
260 // boolean should be set if the parameter is given, etc. may be specified
261 // using the various set* calls of the Xyce::Device::Descriptor class.
262 //
263 // It is important to note that since parameters defined by addPar are
264 // defined as metadata separate from any instance, it is not possible to
265 // establish interrelationships between parameter defaults here. Parameter
266 // defaults specified in addPar are always constants. If a device requires
267 // that parameter defaults depend on values of other parameters, this has to
268 // be done in the instance constructor. Examples of such parameters in this
269 // device arethe "DTEMP" and "W" parameters, which are set to special defaults
270 // at device instantiation time. Defaults for those parameters in the addPar
271 // calls (and hence in the LaTeX tables in the reference guide produced from
272 // this data) are misleading.
273 //
275 {
276  p.addPar("IVRELATION", 0, &MemristorTEAM::Instance::IVrelation)
277  .setUnit(U_NONE)
278  .setDescription("IV relationship to use, 0 is linear, 1 is nonlinear");
279 
280 }
281 
282 //-----------------------------------------------------------------------------
283 // Function : Xyce::Device::MemristorTEAM::Traits::loadModelParameters
284 // Purpose :
285 // Special Notes : The addPar calls here were refactored and moved here
286 // from the model constructor. Those addPars had been
287 // in place from 2/4/2005.
288 // Scope : private
289 // Creator : Richard Schiek, Electrical Models & Simulations
290 // Creation Date : 10/23/2014
291 //-----------------------------------------------------------------------------
292 //
293 // Loads the parameter definition into the model parameter map.
294 //
295 // @param p model parameter map
296 //
297 // @see Xyce::Device::MemristorTEAM::Traits::loadInstanceParameters
298 //
299 //
301 {
302 
303  // Create parameter definitions for parameter member variables
304  p.addPar("KON", -8.0e-13, &MemristorTEAM::Model::kOn_)
305  .setUnit(U_MSM1)
306  .setDescription("Modeling Coefficient");
307  p.addPar("KOFF", 8.0e-13, &MemristorTEAM::Model::kOff_)
308  .setUnit(U_MSM1)
309  .setDescription("Modeling Coefficient");
310  p.addPar("ALPHAON", 3.0, &MemristorTEAM::Model::alphaOn_)
311  .setUnit(U_NONE)
312  .setDescription("Modeling Coefficient");
313  p.addPar("ALPHAOFF", 3.0, &MemristorTEAM::Model::alphaOff_)
314  .setUnit(U_NONE)
315  .setDescription("Modeling Coefficient");
316  p.addPar("XON", 0.0, &MemristorTEAM::Model::xOn_)
317  .setUnit(U_METER)
318  .setDescription("Modeling Coefficient");
319  p.addPar("XOFF", 3.0e-9, &MemristorTEAM::Model::xOff_)
320  .setUnit(U_METER)
321  .setDescription("Modeling Coefficient");
322  p.addPar("RON", 50.0, &MemristorTEAM::Model::ROn_)
323  .setUnit(U_OHM)
324  .setDescription("Resistence in on state");
325  p.addPar("ROFF", 1.0e3, &MemristorTEAM::Model::ROff_)
326  .setUnit(U_OHM)
327  .setDescription("Resistence in off state");
328  p.addPar("ION", 8.9e-6, &MemristorTEAM::Model::iOn_)
329  .setUnit(U_AMP)
330  .setDescription("Current scale in On state");
331  p.addPar("IOFF", 1.15e-4, &MemristorTEAM::Model::iOff_)
332  .setUnit(U_OHM)
333  .setDescription("Current scale in off state");
334  p.addPar("XSCALING", 1.0, &MemristorTEAM::Model::xScaling_)
335  .setUnit(U_NONE)
336  .setDescription("Scaling for x variable. For example 1e9 if x will be in units of nanometers.");
337 
338  p.addPar("D", 1.15e-4, &MemristorTEAM::Model::D_)
339  .setUnit(U_NONE)
340  .setDescription("Window Function Parameter (windows 1, 2 and 3)");
341  p.addPar("P", 1.15e-4, &MemristorTEAM::Model::p_)
342  .setUnit(U_NONE)
343  .setDescription("Window Function Parameter (windows 1, 2 and 3)");
344  p.addPar("J", 1.15e-4, &MemristorTEAM::Model::j_)
345  .setUnit(U_NONE)
346  .setDescription("Window Function Parameter (window 3)");
347 
348  p.addPar("AON", 0.0, &MemristorTEAM::Model::aOn_)
349  .setUnit(U_METER)
350  .setDescription("Window Function Parameter (window 4)");
351  p.addPar("AOFF", 3.0e-9, &MemristorTEAM::Model::aOff_)
352  .setUnit(U_METER)
353  .setDescription("Window Function Parameter (window 4)");
354 
355  p.addPar("WC", 1.07e-12, &MemristorTEAM::Model::wc_)
356  .setUnit(U_METER)
357  .setDescription("Window Function Parameter (window 4)");
358 
360  .setUnit(U_NONE)
361  .setDescription("Type of windowing function: 0-None, 1-Jogelkar, 2-Biolek, 3-Prodromakis, 4-Kvatinsky");
362 
363  p.addPar("RESNOISE", false,&MemristorTEAM::Model::randomResNoiseOn_)
364  .setUnit(U_NONE)
365  .setDescription("RTN model in resistance (on/off)" );
366 
368  .setUnit(U_NONE)
369  .setDescription("RTN model in resistance: seed" );
370 
372  .setUnit(U_NONE)
373  .setDescription("RTN model: lambda" );
374 
376  .setUnit(U_SECOND)
377  .setDescription("RTN model in resistance: Update time" );
378 
380  .setUnit(U_SECOND)
381  .setDescription("RTN model in resistance: Minimum allowed update time" );
382 
384  .setUnit(U_OHM)
385  .setDescription("RTN model in resistance: Base change in resistance for RTN" );
386 
387  p.addPar("RESDELTAGRAD", 0.0,&MemristorTEAM::Model::randomResDeltaGrad_)
388  .setUnit(U_NONE)
389  .setDescription("RTN model in resistance: Base change in resistance for RTN scaled by R" );
390 
391 }
392 
393 //-----------------------------------------------------------------------------
394 // Function : Xyce::Device::MemristorTEAM::Instance::Instance
395 // Purpose : Instance constructor
396 // Special Notes :
397 // Scope : public
398 // Creator : Richard Schiek, Electrical Models & Simulations
399 // Creation Date : 10/23/2014
400 //-----------------------------------------------------------------------------
401 //
402 // Construct a memristor instance.
403 //
404 
406  const Configuration & configuration,
407  const InstanceBlock & instance_block,
408  Model & model,
409  const FactoryBlock & factory_block)
410  : DeviceInstance(instance_block, configuration.getInstanceParameters(), factory_block),
411  model_(model),
412  IVrelation(0),
413  G(0.0),
414  i0(0.0),
415  resNoiseLastUpdateStep(0),
416  resNoiseLastUpdateTime(0.0),
417  resNoiseNextUpdateTime(0.0),
418  resNoiseHiLoState(0),
419  resNoiseRfactor(1.0),
420  li_Pos(-1),
421  li_Neg(-1),
422  li_x(-1),
423  li_store_dev_i(-1),
424  li_store_R(-1),
425  li_store_tdt(-1),
426  li_branch_data(0),
427  APosEquPosNodeOffset(-1),
428  APosEquNegNodeOffset(-1),
429  APosEquXNodeOffset(-1),
430  ANegEquPosNodeOffset(-1),
431  ANegEquNegNodeOffset(-1),
432  ANegEquXNodeOffset(-1),
433  XEquVPosOffset(-1),
434  XEquVNegOffset(-1),
435  XEquXOffset(-1)
437  ,
438  f_PosEquPosNodePtr(0),
439  f_PosEquNegNodePtr(0),
440  f_PosEquXNodePtr(0),
441  f_NegEquPosNodePtr(0),
442  f_NegEquNegNodePtr(0),
443  f_NegEquXNodePtr(0),
444  f_XEquPosNodePtr(0),
445  f_XEquNegNodePtr(0),
446  f_XEquXNodePtr(0),
447  q_XEquXNodePtr(0)
448  // ,
449 #endif
450  //resSens(*this)
451 {
452  // Initialize DeviceInstance values
453  numIntVars = 1; // Initialize number if internal nodes in DeviceInstance
454  numExtVars = 2; // Initialize number if external nodes in DeviceInstance
455  numStateVars = 0; // Initialize number if state variables in DeviceInstance
456  setNumStoreVars(2); // Initialize number if store variables in DeviceInstance
457  numLeadCurrentStoreVars = 1; // Initialize number if lead current variables
458  //in DeviceInstance
459 
460  setNumBranchDataVars(0); // by default don't allocate space in branch vectors
461  numBranchDataVarsIfAllocated = 1; // this is the space to allocate if lead current or power is needed.
462 
464 
465  // Set params to constant default values from parameter definition
467 
468  // Set params according to instance line and constant defaults from metadata
469  setParams(instance_block.params);
470 
471  // Calculate any parameters specified as expressions
473 
474  // Process the parameters to complete initialization
475  processParams();
476 
477  // if the random resistance noise model is on the update the "update" time.
479  {
480  // get a poisson distributed number with mean lambda
481  //int aNum = model_.randomNumberGen_->poissonRandom(model_.randomResNoiseLambda_);
482  // the actual update interval is the poisson distributed number times the base update time.
483  //resNoiseNextUpdateTime = aNum * model_.randomResUpdateTime_ + model_.randomResEpsilonUpdateTime_;
485  // don't allow a zero time for the next update.
486  //if( aNum == 0 )
487  // resNoiseNextUpdateTime = model_.randomResEpsilonUpdateTime_;
488  }
489 
490 
491 }
492 
493 //-----------------------------------------------------------------------------
494 // Function : Xyce::Device::MemristorTEAM::Instance::processParams
495 // Purpose :
496 // Special Notes :
497 // Scope : public
498 // Creator : Richard Schiek, Electrical Models & Simulations
499 // Creation Date : 10/23/2014
500 //-----------------------------------------------------------------------------
501 // Process parameters.
502 //
503 // @return true on success
504 //
505 // In general, the processParams method is intended as a place for complex
506 // manipulation of parameters that must happen if temperature or other
507 // external changes happen.
508 //
510 {
511  // now set the temperature related stuff.
512  double temp=0;
513  return updateTemperature(temp);
514 }
515 
516 //-----------------------------------------------------------------------------
517 // Function : Xyce::Device::MemristorTEAM::Instance::registerLIDs
518 // Purpose :
519 // Special Notes :
520 // Scope : public
521 // Creator : Richard Schiek, Electrical Models & Simulations
522 // Creation Date : 10/23/2014
523 //-----------------------------------------------------------------------------
524 //
525 // Register local IDs
526 //
527 // Register the local internal and external node IDs.
528 //
529 // @param intLIDVecRef internal local IDs from topology package
530 // @param extLIDVecRef external local IDs from topology package
531 //
533  const std::vector<int> & intLIDVecRef,
534  const std::vector<int> & extLIDVecRef)
535 {
536  AssertLIDs(intLIDVecRef.size() == numIntVars);
537  AssertLIDs(extLIDVecRef.size() == numExtVars);
538 
539  // Copy the local ID lists.
540  intLIDVec = intLIDVecRef; // Set the internal local IDs in DeviceInstance
541  extLIDVec = extLIDVecRef; // Set the external local IDs in DeviceInstance
542 
543  li_Pos = extLIDVec[0];
544  li_Neg = extLIDVec[1];
545 
546  // This lid is for the internal variable for layer thickess that determines resistence
547  li_x = intLIDVec[0];
548 
549 }
550 
551 //-----------------------------------------------------------------------------
552 // Function : Xyce::Device::MemristorTEAM::Instance::registerStateLIDs
553 // Purpose :
554 // Special Notes :
555 // Scope : public
556 // Creator : Richard Schiek, Electrical Models & Simulations
557 // Creation Date : 10/23/2014
558 //-----------------------------------------------------------------------------
559 //
560 // Register the local state IDs
561 //
562 // @note The memristor does not have any state vars, so this function
563 // does nothing.
564 //
565 // @param staLIDVecRef State variable local IDs
566 //
567 // In general, devices may declare at construction that they require storage
568 // locations in the "state vector." Topology assigns locations in the
569 // state vector and returns "Local IDs" (LIDs) for devices to use for their
570 // state vector entries. If a device uses state vector entries, it
571 // uses the registerStateLIDs method to save the local IDs for later use.
572 //
573 // @author Robert Hoekstra, SNL, Parallel Computational Sciences
574 // @date 06/12/02
575 
576 void Instance::registerStateLIDs(const std::vector<int> & staLIDVecRef)
577 {
578  AssertLIDs(staLIDVecRef.size() == numStateVars);
579 }
580 
581 //-----------------------------------------------------------------------------
582 // Function : Xyce::Device::MemristorTEAM::Instance::registerStoreLIDs
583 // Purpose :
584 // Special Notes :
585 // Scope : public
586 // Creator : Richard Schiek, Electrical Models & Simulations
587 // Creation Date : 12/18/2012
588 //-----------------------------------------------------------------------------
589 // Register the local store IDs
590 //
591 // In addition to state vector, Xyce maintains a separate datastructure
592 // called a "store" vector. As with other such vectors, the device
593 // declares at construction time how many store vector entries it needs,
594 // and later Topology assigns locations for devices, returning LIDs.
595 //
596 //
597 //
598 
599 void Instance::registerStoreLIDs(const std::vector<int> & stoLIDVecRef)
600 {
601  AssertLIDs(stoLIDVecRef.size() == getNumStoreVars());
602  li_store_R = stoLIDVecRef[0];
603  li_store_tdt = stoLIDVecRef[1];
604  if (loadLeadCurrent)
605  {
606  li_store_dev_i = stoLIDVecRef[2];
607  }
608 
609 }
610 
611 //-----------------------------------------------------------------------------
612 // Function : Xyce::Device::MemristorTEAM::Instance::registerBranchDataLIDs
613 // Purpose :
614 // Special Notes :
615 // Scope : public
616 // Creator : Richard Schiek, Electrical Systems Modeling
617 // Creation Date : 12/18/2012
618 //-----------------------------------------------------------------------------
619 // Register the local store IDs
620 //
621 // In addition to state vector, Xyce maintains a separate datastructure
622 // called a "branch data" vector. As with other such vectors, the device
623 // declares at construction time how many branch vector entries it needs,
624 // and later Topology assigns locations for devices, returning LIDs.
625 //
626 // These LIDs are stored in this method for later use.
627 //
628 // The memristor device uses exactly one "branch data vector" element, where
629 // it keeps the "lead current" that may be used on .PRINT lines as
630 // "I(ymemristor)" for the current through ymemristor. and a junction voltage.
631 //
632 // @param stoLIDVecRef Store variable local IDs
633 //
634 // @author Richard Schiek, Electrical Systems Modeling
635 // @date 12/18/2012
636 
637 void Instance::registerBranchDataLIDs(const std::vector<int> & branchLIDVecRef)
638 {
639  AssertLIDs(branchLIDVecRef.size() == getNumBranchDataVars());
640 
641  if (loadLeadCurrent)
642  {
643  li_branch_data= branchLIDVecRef[0];
644  }
645 }
646 
647 
648 
649 //-----------------------------------------------------------------------------
650 // Function : Instance::loadNodeSymbols
651 // Purpose :
652 // Special Notes :
653 // Scope : public
654 // Creator : Richard Schiek, Electrical Sysetms Modeling
655 // Creation Date : 05/13/05
656 //-----------------------------------------------------------------------------
657 void Instance::loadNodeSymbols(Util::SymbolTable &symbol_table) const
658 {
659  addInternalNode(symbol_table, li_x, getName(), "x");
660 
661  addStoreNode(symbol_table, li_store_R, getName(), "R");
662  addStoreNode(symbol_table, li_store_tdt, getName(), "TDT");
663 
664  if (loadLeadCurrent)
665  {
666  addBranchDataNode( symbol_table, li_branch_data, getName(), "BRANCH_D");
667  }
668 }
669 
670 
671 
672 
673 //-----------------------------------------------------------------------------
674 // Function : Xyce::Device::MemristorTEAM::Instance::registerJacLIDs
675 // Purpose :
676 // Special Notes :
677 // Scope : public
678 // Creator : Richard Schiek, Electrical Models & Simulations
679 // Creation Date : 10/23/2014
680 //-----------------------------------------------------------------------------
681 //
682 // Register the Jacobian local IDs
683 //
684 // @param jacLIDVec Jacobian local Ids
685 //
686 // @see Xyce::Device::MemristorTEAM::Instance::initializeJacobianStamp
687 //
688 // Having established local IDs for the solution variables, Topology must
689 // also assign local IDs for the elements of the Jacobian matrix.
690 //
691 // For each non-zero element that was identified in the jacobianStamp,
692 // Topology will assign a Jacobian LID. The jacLIDVec will have the
693 // same structure as the JacStamp, but the values in the jacLIDVec will
694 // be offsets into the row of the sparse Jacobian matrix corresponding
695 // to the non-zero identified in the stamp.
696 //
697 // These offsets are stored in this method for use later when we load
698 // the Jacobian.
699 //
700 // @author Robert Hoekstra, SNL, Parallel Computational Sciences
701 // @date 08/27/01
702 
703 void Instance::registerJacLIDs(const std::vector< std::vector<int> > & jacLIDVec)
704 {
705  // Let DeviceInstance do its work.
707 
708  // Store offsets of the components of the Jacobian of this instance
709  APosEquPosNodeOffset = jacLIDVec[0][0];
710  APosEquNegNodeOffset = jacLIDVec[0][1];
711  APosEquXNodeOffset = jacLIDVec[0][2];
712  ANegEquPosNodeOffset = jacLIDVec[1][0];
713  ANegEquNegNodeOffset = jacLIDVec[1][1];
714  ANegEquXNodeOffset = jacLIDVec[1][2];
715  XEquVPosOffset = jacLIDVec[2][0];
716  XEquVNegOffset = jacLIDVec[2][1];
717  XEquXOffset = jacLIDVec[2][2];
718 
719 }
720 
721 //-----------------------------------------------------------------------------
722 // Function : Xyce::Device::MemristorTEAM::Instance::setupPointers
723 // Purpose :
724 // Special Notes :
725 // Scope : public
726 // Creator : Richard Schiek, Electrical Models & Simulations
727 // Creation Date : 10/23/2014
728 //-----------------------------------------------------------------------------
729 //
730 // Setup direct access pointer to solution matrix and vectors.
731 //
732 // @see Xyce::Device::MemristorTEAM::Instance::registerJacLIDs
733 //
734 // As an alternative to the row offsets defined in registerJacLIDs, it
735 // is also possible to obtain direct pointers of the Jacobian elements.
736 //
737 // This method uses the offsets obtained in registerJacLIDs to retrieve
738 // the pointers.
739 //
740 // In this device the pointers to the matrix are only saved
741 // (and are only used for matrix access) if
742 // Xyce_NONPOINTER_MATRIX_LOAD is NOT defined at compile time. For
743 // some devices the use of pointers instead of array indexing can be
744 // a performance enhancement.
745 //
746 // Use of pointers in this device is disabled by defining
747 // Xyce_NONPOINTER_MATRIX_LOAD at compile time. When disabled, array
748 // indexing with the offsets from registerJacLIDs is used in
749 // the matrix load methods.
750 //
751 // @author Eric Keiter, SNL
752 // @date 11/30/08
753 
755 {
756 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
757  Linear::Matrix & dFdx = *(extData.dFdxMatrixPtr);
758  Linear::Matrix & dQdx = *(extData.dQdxMatrixPtr);
767  f_XEquXNodePtr = &(dFdx[li_x][XEquXOffset]);
768  q_XEquXNodePtr = &(dQdx[li_x][XEquXOffset]);
769 #endif
770 }
771 
772 //-----------------------------------------------------------------------------
773 // Function : Xyce::Device::MemristorTEAM::Instance::updateIntermediateVars
774 // Purpose :
775 // Special Notes :
776 // Scope : public
777 // Creator : Richard Schiek, Electrical Models & Simulations
778 // Creation Date : 10/23/2014
779 //-----------------------------------------------------------------------------
780 //
781 // Update the intermediate variables
782 //
783 // @return true on success
784 //
785 // The bulk of any device's computation is handled in the instance class'
786 // updateIntermediateVars method. For the memristor, this is
787 // merely the computation of the current through the device given the
788 // voltage difference between its terminals.
789 //
790 // Intermediate variables computed here are used in the methods that
791 // load data into the F, Q, dFdX and dQdX data structures.
792 //
793 // @note This method is called by the updatePrimaryState
794 // method. Since the MemristorTEAM class reimplements the "Master"
795 // "loadState" function that loads the contributions from all
796 // such devices in a single loop, THIS FUNCTION IS NOT ACTUALLY
797 // USED!
798 //
799 //
801 {
802  double * solVec = extData.nextSolVectorRawPtr;
803  double v_pos = solVec[li_Pos];
804  double v_neg = solVec[li_Neg];
805  double x = solVec[li_x];
806  double Reff = model_.ROn_;
807  //
808  // Calclate Reff and deriative term dReff/dx
809  //
810  {
811  // extra scope to limit variable extent. Otherwise one will have lots of
812  // similar variables
813  Sacado::Fad::SFad<double,1> varX( 1, 0, x );
814  Sacado::Fad::SFad<double,1> paramRon( model_.ROn_ );
815  Sacado::Fad::SFad<double,1> paramRoff( model_.ROff_ );
816  Sacado::Fad::SFad<double,1> paramXon( model_.xOn_ );
817  Sacado::Fad::SFad<double,1> paramXoff( model_.xOff_ );
818  Sacado::Fad::SFad<double,1> resultFad;
819  if(IVrelation==0)
820  {
821  // liner current voltage relationship
822  ReffLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
823 
824  }
825  else if( IVrelation==1 )
826  {
827  ReffNonLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
828  }
829  Reff = resultFad.val();
830  dReffdx = resultFad.dx(0);
831  dReffdvpos = -v_pos * pow( Reff, 2 ) * dReffdx;
832  dReffdvneg = v_neg * pow( Reff, 2 ) * dReffdx;
833  }
834 
835  // if the random noise in the resistance model is on
836  // and this is the start of a new step, then get a new
837  // gaussian distributed variable
838  double rfactor=1.0;
839  if( model_.randomResNoiseOn_ && (getSolverState().timeStepNumber_ != resNoiseLastUpdateStep) )
840  {
843  }
844 
845  G=1.0/(rfactor*Reff);
846  i0= (v_pos-v_neg)*G;
847 
848  {
849  // calculate the xVar F equation update and the
850  // derivatives d(F(x))/dVpos, d(F(x))/dVneg and d(F(x))/dx
851  Sacado::Fad::SFad<double,3> varVpos( 3, 0, v_pos);
852  Sacado::Fad::SFad<double,3> varVneg( 3, 1, v_neg);
853  Sacado::Fad::SFad<double,3> varX( 3, 2, x );
854  Sacado::Fad::SFad<double,3> parami( i0 );
855  Sacado::Fad::SFad<double,3> paramG( G );
856  Sacado::Fad::SFad<double,3> paramiOff( model_.iOff_ );
857  Sacado::Fad::SFad<double,3> paramiOn( model_.iOn_ );
858  Sacado::Fad::SFad<double,3> paramkOff( model_.kOff_ );
859  Sacado::Fad::SFad<double,3> paramkOn( model_.kOn_ );
860  Sacado::Fad::SFad<double,3> paramAlphaOff( model_.alphaOff_ );
861  Sacado::Fad::SFad<double,3> paramAlphaOn( model_.alphaOn_ );
862  Sacado::Fad::SFad<double,3> resultFad;
863 
864  xVarFterm(varVpos, varVneg, varX, paramG, paramiOff, paramiOn, paramkOff, paramkOn,
865  paramAlphaOff, paramAlphaOn, resultFad);
866 
867  // calculate the window function value
868  Sacado::Fad::SFad<double,3> windowFunctionFad;
869 
870  if( model_.windowType_ == 0)
871  {
872  // no window
873  windowFunctionFad = 1.0;
874  }
875  else if( model_.windowType_ == 1)
876  {
877  Sacado::Fad::SFad<double,3> paramD( model_.D_ );
878  Sacado::Fad::SFad<double,3> paramP( model_.p_ );
879  JogelkarWindowFunction( varX, paramD, paramP, windowFunctionFad );
880  }
881  else if( model_.windowType_ == 2)
882  {
883  Sacado::Fad::SFad<double,3> paramD( model_.D_ );
884  Sacado::Fad::SFad<double,3> paramP( model_.p_ );
885  BiolekWindowFunction( varX, paramD, paramP, parami, windowFunctionFad );
886  }
887  else if( model_.windowType_ == 3)
888  {
889  Sacado::Fad::SFad<double,3> paramD( model_.D_ );
890  Sacado::Fad::SFad<double,3> paramP( model_.p_ );
891  Sacado::Fad::SFad<double,3> paramJ( model_.j_ );
892  ProdromakisWindowFunction( varX, paramD, paramP, paramJ, windowFunctionFad );
893  }
894  else if( model_.windowType_ == 4)
895  {
896  Sacado::Fad::SFad<double,3> paramAOn( model_.aOn_ );
897  Sacado::Fad::SFad<double,3> paramAOff( model_.aOff_ );
898  Sacado::Fad::SFad<double,3> paramWc( model_.wc_ );
899  TEAMWindowFunctionF( varX, parami, paramAOff, paramAOn, paramWc, windowFunctionFad );
900  }
901 
902  xVarFContribution = resultFad.val();
903  dxFEqdVpos = resultFad.dx(0);
904  dxFEqdVneg = resultFad.dx(1);
905  dxFEqdx = resultFad.dx(2);
906  //xVarFContribution = 1.0;
907  //dxFEqdVpos = 0.0;
908  //dxFEqdVneg = 0.0;
909  //dxFEqdx = 0.0;
910 
911  }
912 
913  //i0 = (v_pos-v_neg)*G;
914 
915  // get derivative terms for jacobian load
916  //
917  // APosEquPosNodeOffset --> G
918  // APosEquNegNodeOffset --> -G
919  // APosEquXNodeOffset --> -Reff^2 dR/dx
920  // ANegEquPosNodeOffset --> -G
921  // ANegEquNegNodeOffset --> G
922  // ANegEquXNodeOffset --> Reff^2 dR/dx
923  // XEquVPosOffset --> d xVarFterm / dVpos
924  // XEquVNegOffset --> d xVarFterm / dVneg
925  // XEquXOffset --> d xVarFterm / dx
926 
927  return true;
928 }
929 
930 //-----------------------------------------------------------------------------
931 // Function : Xyce::Device::MemristorTEAM::Instance::updatePrimaryState
932 // Purpose :
933 // Special Notes :
934 // Scope : public
935 // Creator : Richard Schiek, Electrical Models & Simulations
936 // Creation Date : 10/23/2014
937 //-----------------------------------------------------------------------------
938 //
939 // Update the state variables.
940 //
941 // @return true on success
942 //
943 // This function is the function that is called from the device manager
944 // each time a load of vectors and Jacobian is required. Its first
945 // job is to call updateIntermediateVars.
946 //
947 // After calling updateIntermediateVars, the updatePrimaryState method
948 // may load state vector elements as needed.
949 //
950 // The memristor device has no state vector elements, so all this method does
951 // in the memristor instance class is call updateIntermediateVars.
952 //
953 // There is no longer a "secondary" state. The "primary" in
954 // this method's name is purely historical.
955 //
956 // @note This method is called by the default implementation of the
957 // loadState master function. Since the MemristorTEAM class reimplements
958 // the "Master" "loadState" function that loads the contributions
959 // from all memristor devices in a single loop, THIS FUNCTION IS NOT
960 // ACTUALLY USED, NOR is the updateIntermediateVars method it calls!
961 // The updatePrimaryState method is only called when a device class
962 // does not re-implement the master class. This can be a source of
963 // confusion when attempting to modify the MemristorTEAM device, or any
964 // other device That reimplements the Master classes.
965 //
966 // @see Xyce::Device::MemristorTEAM::Master::updateState
967 //
968 // @author Eric Keiter, SNL, Parallel Computational Sciences
969 // @date 01/29/01
970 //
972 {
973  return updateIntermediateVars();
974 }
975 
977 {
978  double * qVec = extData.daeQVectorRawPtr;
979  double * solVec = extData.nextSolVectorRawPtr;
980  double x = solVec[li_x];
981 
982  qVec[li_x] += x;
983 
984  return true;
985 }
986 
987 
988 //-----------------------------------------------------------------------------
989 // Function : Xyce::Device::MemristorTEAM::Instance::loadDAEFVector
990 // Purpose :
991 // Special Notes :
992 // Scope : public
993 // Creator : Richard Schiek, Electrical Models & Simulations
994 // Creation Date : 10/23/2014
995 //-----------------------------------------------------------------------------
996 //
997 // Load the DAE F vector.
998 //
999 // @return true on success
1000 //
1001 // The Xyce DAE formulation solves the differential-algebraic
1002 // equations \f$F(x)+dQ(x)/dt=0\f$ These are vector equations
1003 // resulting from the modified nodal analysis of the circuit.
1004 //
1005 // This method loads the F-vector contributions for a single memristor
1006 // instance.
1007 //
1008 // In this method, the offsets defined in registerLIDs are used to
1009 // store the device's F contributions into the F vector.
1010 //
1011 // The Q vector is used for devices that store charge or magnetic
1012 // energy, and since the memristor does none of that, the F vector
1013 // contributions are the whole of the memristor's contribution to the
1014 // full set of DAEs.
1015 //
1016 // @note This method is called by the default implementation of the
1017 // loadDAEVectors master function. Since the MemristorTEAM class
1018 // reimplements the "Master" "loadDAEVectors" function that loads the
1019 // contributions from all memristor devices in a single loop, THIS
1020 // FUNCTION IS NOT ACTUALLY USED. The loadDAEFVector method is only
1021 // called when a device class does not re-implement the master class.
1022 // This can be a source of confusion when attempting to modify the MemristorTEAM
1023 // device, or any other device that reimplements the Master classes.
1024 //
1025 // @see Xyce::Device::MemristorTEAM::Master::loadDAEVectors
1026 //
1027 // @author Eric Keiter, SNL, Parallel Computational Sciences
1028 // @date 01/24/03
1029 //
1031 {
1032  double * fVec = extData.daeFVectorRawPtr;
1033  fVec[li_Pos] += i0;
1034  fVec[li_Neg] += -i0;
1035  fVec[li_x] += xVarFContribution;
1036 
1037  if( loadLeadCurrent )
1038  {
1039  double * leadF = extData.nextLeadCurrFCompRawPtr;
1040  double * junctionV = extData.nextJunctionVCompRawPtr;
1041  double * solVec = extData.nextSolVectorRawPtr;
1042  leadF[li_branch_data] = i0;
1043  junctionV[li_branch_data] = solVec[li_Pos] - solVec[li_Neg];
1044  }
1045 
1046 
1047  return true;
1048 }
1049 
1051 {
1052  Linear::Matrix & dQdx = *(extData.dQdxMatrixPtr);
1053  dQdx[li_x][XEquXOffset] = 1.0;
1054  return true;
1055 }
1056 
1057 
1058 //-----------------------------------------------------------------------------
1059 // Function : Xyce::Device::MemristorTEAM::Instance::loadDAEdFdx
1060 // Purpose : Loads the F-vector contributions for a single
1061 // memristor instance.
1062 // Special Notes :
1063 // Scope : public
1064 // Creator : Richard Schiek, Electrical Models & Simulations
1065 // Creation Date : 10/23/2014
1066 //-----------------------------------------------------------------------------
1067 //
1068 // Load the DAE the derivative of the F vector with respect to the
1069 // solution vector x, dFdx
1070 //
1071 // Loads the contributions for a single memristor instance to the
1072 // dFdx matrix (the F contribution to the Jacobian).
1073 //
1074 // This method uses the Jacobian LIDs (row offsets) that were stored by
1075 // registerJacLIDs.
1076 //
1077 // @see Xyce::Device::MemristorTEAM::Instance::registerJacLIDs
1078 //
1079 // @note This method is called by the default implementation of the
1080 // loadDAEMatrices master function. Since the MemristorTEAM class
1081 // reimplements the "Master" "loadDAEMatrices" function that loads the
1082 // contributions from all memristor devices in a single loop, THIS
1083 // FUNCTION IS NOT ACTUALLY USED. The loadDAEdFdx method is only
1084 // called when a device class does not re-implement the master class.
1085 // This can be a source of confusion when attempting to modify the MemristorTEAM
1086 // device, or any other device that reimplements the Master classes.
1087 //
1088 // @see Xyce::Device::MemristorTEAM::Master::loadDAEMatrices
1089 //
1090 // @return true on success
1091 //
1092 // @author Eric Keiter, SNL, Parallel Computational Sciences
1093 // @date 03/05/04
1095 {
1096  Linear::Matrix & dFdx = *(extData.dFdxMatrixPtr);
1097  dFdx[li_Pos][APosEquPosNodeOffset] += G;
1098  dFdx[li_Pos][APosEquNegNodeOffset] -= G;
1100  dFdx[li_Neg][ANegEquPosNodeOffset] -= G;
1101  dFdx[li_Neg][ANegEquNegNodeOffset] += G;
1103  dFdx[li_x][XEquVPosOffset] += dxFEqdVpos;
1104  dFdx[li_x][XEquVNegOffset] += dxFEqdVneg;
1105  dFdx[li_x][XEquXOffset] += dxFEqdx;
1106  return true;
1107 }
1108 
1109 //-----------------------------------------------------------------------------
1110 // Function : Xyce::Device::MemristorTEAM::Instance::updateTemperature
1111 // Purpose :
1112 // Special Notes :
1113 // Scope : public
1114 // Creator : Richard Schiek, Electrical Models & Simulations
1115 // Creation Date : 10/23/2014
1116 //-----------------------------------------------------------------------------
1117 //
1118 // Update the parameters that depend on the temperature of the device
1119 //
1120 // @param temp_tmp temperature
1121 //
1122 // Xyce has a number of mechanisms that allow temperature to be changed
1123 // after a device has been instantiated. These include .STEP loops over
1124 // temperature. When temperature is changed, any device that has parameters
1125 // that depend on temperature must be updated. That updating happens here.
1126 //
1127 // The MemristorTEAM device supports temperature-dependent resistance through its
1128 // TC1 (linear dependence) and TC2 (quadratic dependence) parameters.
1129 // If these parameters are specified, the resistance must be updated.
1130 //
1131 // @return true on success
1132 //
1133 // @author Tom Russo, Component Information and Models
1134 // @date 02/27/01
1135 bool Instance::updateTemperature(const double & temp_tmp)
1136 {
1137  bool bsuccess = true;
1138  /*
1139  double difference, factor;
1140 
1141  if (temp_tmp != -999.0)
1142  temp = temp_tmp;
1143  difference = temp - model_.tnom;
1144  factor = model_.resistanceMultiplier*(1.0 + tempCoeff1*difference + tempCoeff2*difference*difference);
1145 
1146  if (R*factor != 0.0)
1147  G = 1.0/(R * factor);
1148  else
1149  G = 0.0;
1150  */
1151  return bsuccess;
1152 }
1153 
1154 //-----------------------------------------------------------------------------
1155 // Function : Xyce::Device::MemristorTEAM::Model::processParams
1156 // Purpose :
1157 // Special Notes :
1158 // Scope : public
1159 // Creator : Richard Schiek, Electrical Models & Simulations
1160 // Creation Date : 10/23/2014
1161 //-----------------------------------------------------------------------------
1162 //
1163 // Process model parameters
1164 //
1165 // @return true on success
1166 //
1167 // @author Eric Keiter, SNL, Parallel Computational Sciences
1168 // @date 6/03/02
1170 {
1171  return true;
1172 }
1173 
1174 //----------------------------------------------------------------------------
1175 // Function : Xyce::Device::MemristorTEAM::Model::processInstanceParams
1176 // Purpose :
1177 // Special Notes :
1178 // Scope : public
1179 // Creator : Richard Schiek, Electrical Models & Simulations
1180 // Creation Date : 10/23/2014
1181 //----------------------------------------------------------------------------
1182 //
1183 // Process the instance parameters of instance owned by this model
1184 //
1185 // This method simply loops over all instances associated with this
1186 // model and calls their processParams method.
1187 //
1188 // @return true
1189 //
1190 // @author Dave Shirely, PSSI
1191 // @date 03/23/06
1192 
1194 {
1195  for (InstanceVector::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
1196  {
1197  (*it)->processParams();
1198  }
1199 
1200  return true;
1201 }
1202 //-----------------------------------------------------------------------------
1203 // Function : Xyce::Device::MemristorTEAM::Model
1204 // Purpose : model block constructor
1205 // Special Notes :
1206 // Scope : public
1207 // Creator : Richard Schiek, Electrical Models & Simulations
1208 // Creation Date : 10/23/2014
1209 //-----------------------------------------------------------------------------
1210 //
1211 // Construct a memristor model from a "model block" that was created
1212 // by the netlist parser.
1213 //
1214 // @param configuration
1215 // @param model_block
1216 // @param factory_block
1217 //
1218 // @author Eric Keiter, SNL, Parallel Computational Sciences
1219 // @date 5/16/00
1221  const Configuration & configuration,
1222  const ModelBlock & model_block,
1223  const FactoryBlock & factory_block)
1224  : DeviceModel(model_block, configuration.getModelParameters(), factory_block),
1225  kOn_(0.0),
1226  kOff_(0.0),
1227  alphaOn_(0.0),
1228  alphaOff_(0.0),
1229  xOn_(0.0),
1230  xOff_(0.0),
1231  randomResNoiseOn_(false),
1232  randomResNoiseSeed_(0),
1233  randomResNoiseLambda_(0.0),
1234  randomResNoiseMean_(0.0),
1235  randomResNoiseSD_(0.0),
1236  randomResUpdateTime_(0.0),
1237  randomResEpsilonUpdateTime_(0.0),
1238  randomResDelta_(0.0),
1239  randomResDeltaGrad_(0.0)
1240 
1241 {
1242  // Set params to constant default values.
1243  setDefaultParams();
1244 
1245  // Set params according to .model line and constant defaults from metadata.
1246  setModParams(model_block.params);
1247 
1248  // Set any non-constant parameter defaults.
1249  //if (!given("TNOM"))
1250  // tnom = getDeviceOptions().tnom;
1251 
1252  // Calculate any parameters specified as expressions.
1254 
1255  // calculate dependent (ie computed) params and check for errors.
1256  processParams();
1257 
1258  // if the random noise in the resistance model is onj
1259  // seed the random number generator
1260  if( randomResNoiseOn_ )
1261  {
1262  randomNumberGen_ = new Xyce::Util::RandomNumbers( randomResNoiseSeed_ );
1263  }
1264 
1265 }
1266 
1267 //-----------------------------------------------------------------------------
1268 // Function : Xyce::Device::MemristorTEAM::Model::~Model
1269 // Purpose : destructor
1270 // Special Notes :
1271 // Scope : public
1272 // Creator : Richard Schiek, Electrical Models & Simulations
1273 // Creation Date : 10/23/2014
1274 //-----------------------------------------------------------------------------
1275 //
1276 // Destroy this model.
1277 //
1278 // Also destroys all instances that use this model.
1279 //
1280 // @author Eric Keiter, SNL, Parallel Computational Sciences
1281 // @date 3/16/00
1283 {
1284  // Destory all owned instances
1285  for (InstanceVector::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
1286  {
1287  delete (*it);
1288  }
1289  // if the random noise in the resistance model is onj
1290  // delete the random number generator
1291  if( randomResNoiseOn_ )
1292  {
1293  delete randomNumberGen_;
1294  }
1295 }
1296 
1297 //-----------------------------------------------------------------------------
1298 // Function : N_DEV_MemristorTEAMModel::printOutInstances
1299 // Purpose : debugging tool.
1300 // Special Notes :
1301 // Scope : public
1302 // Creator : Richard Schiek, Electrical Models & Simulations
1303 // Creation Date : 10/23/2014
1304 //-----------------------------------------------------------------------------
1305 //
1306 // Print instances associated with this model.
1307 //
1308 // Used only for debugging
1309 //
1310 // @param os output stream
1311 //
1312 // @return reference to output stream
1313 //
1314 // @author Eric Keiter, SNL, Parallel Computational Sciences
1315 // @date 4/03/00
1316 //
1317 std::ostream &Model::printOutInstances(std::ostream &os) const
1318 {
1319  os << std::endl;
1320  os << "Number of MemristorTEAM Instances: " << instanceContainer.size() << std::endl;
1321  os << " name model name Parameters" << std::endl;
1322 
1323  int i = 0;
1324  for (InstanceVector::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
1325  {
1326  os << " " << i << ": " << (*it)->getName() << "\t";
1327  os << getName();
1328  //os << "\t\tR(Tnom) = " << (*it)->R;
1329  os << "\tG(T) = " << (*it)->G;
1330  os << std::endl;
1331  ++i;
1332  }
1333 
1334  os << std::endl;
1335 
1336  return os;
1337 }
1338 
1339 //-----------------------------------------------------------------------------
1340 // Function : Model::forEachInstance
1341 // Purpose :
1342 // Special Notes :
1343 // Scope : public
1344 // Creator : Richard Schiek, Electrical Models & Simulations
1345 // Creation Date : 10/23/2014
1346 //-----------------------------------------------------------------------------
1347 // Apply a device instance "op" to all instances associated with this
1348 // model
1349 //
1350 // @param[in] op Operator to apply to all instances.
1351 //
1352 //
1353 void Model::forEachInstance(DeviceInstanceOp &op) const /* override */
1354 {
1355  for (std::vector<Instance *>::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
1356  op(*it);
1357 }
1358 
1359 //-----------------------------------------------------------------------------
1360 // Function : Xyce::Device::MemristorTEAM::Master::updateState
1361 // Purpose :
1362 // Special Notes :
1363 // Scope : public
1364 // Creator : Richard Schiek, Electrical Models & Simulations
1365 // Creation Date : 10/23/2014
1366 //-----------------------------------------------------------------------------
1367 //
1368 // Update state for all memristor instances, regardless of model.
1369 //
1370 // @param solVec solution vector
1371 // @param staVec state vector
1372 // @param stoVec store vector
1373 //
1374 // @return true on success
1375 //
1376 // @note Because the memristor device re-implements the base-class
1377 // Master::updateState, the Instance::updatePrimaryState method is never
1378 // called, nor is the Instance::updateIntermediateVars method. This method
1379 // replaces those, and does the same work but inside a loop over all
1380 // memristor instances.
1381 //
1382 // @see Xyce::Device::MemristorTEAM::Instance::updatePrimaryState
1383 // @author Eric Keiter, SNL
1384 // @date 11/26/08
1385 
1386 bool Master::updateState(double * solVec, double * staVec, double * stoVec)
1387 {
1388  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1389  {
1390  Instance & ri = *(*it);
1391 
1392  double v_pos = solVec[ri.li_Pos];
1393  double v_neg = solVec[ri.li_Neg];
1394  double x = solVec[ri.li_x];
1395  double Reff = ri.model_.ROn_;
1396  //
1397  // Calclate Reff and deriative term dReff/dx
1398  //
1399  {
1400  // extra scope to limit variable extent. Otherwise one will have lots of
1401  // similar variables
1402  Sacado::Fad::SFad<double,1> varX( 1, 0, x );
1403  Sacado::Fad::SFad<double,1> paramRon( ri.resNoiseRfactor*ri.model_.ROn_ );
1404  Sacado::Fad::SFad<double,1> paramRoff( ri.resNoiseRfactor*ri.model_.ROff_ );
1405  Sacado::Fad::SFad<double,1> paramXon( ri.model_.xOn_ * ri.model_.xScaling_ );
1406  Sacado::Fad::SFad<double,1> paramXoff( ri.model_.xOff_ * ri.model_.xScaling_ );
1407  Sacado::Fad::SFad<double,1> resultFad;
1408  if(ri.IVrelation==0)
1409  {
1410  // liner current voltage relationship
1411  //Xyce::dout() << ri.getName() << "About to call ReffLin..." << ri.model_.ROn_ << " , " << ri.model_.ROff_ << ", " << ri.model_.xOn_ << ", " << ri.model_.xOff_ << std::endl;
1412 
1413  ReffLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
1414  //fval = Ron + (X - Xon) * (Roff-Ron)/(Xoff - Xon);
1415  //dReff/dV1 = 0
1416  //dReff/dV2 = 0;
1417  //dReff/dx = (Roff-Ron)/(Xoff - Xon);
1418  }
1419  else if( ri.IVrelation==1 )
1420  {
1421  //Xyce::dout() << ri.getName() << "About to call ReffNonLin..." << std::endl;
1422  ReffNonLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
1423  }
1424  ri.Reff = resultFad.val();
1425  ri.dReffdx = resultFad.dx(0);
1426  ri.dReffdvpos = (v_pos-v_neg)*(-pow( ri.Reff, 2 ) * ri.dReffdx) / (ri.model_.xScaling_);
1427  ri.dReffdvneg = -(v_pos-v_neg)*(-pow( ri.Reff, 2 ) * ri.dReffdx) / (ri.model_.xScaling_);
1428  //ri.dReffdvpos = 0.0;
1429  //ri.dReffdvneg = 0.0;
1430 
1431  }
1432 
1433  ri.G = 1.0/(ri.Reff);
1434  ri.i0 = (v_pos-v_neg)*ri.G;
1435 
1436  // if the random noise in the resistance model is on, then
1437  // with the estimate for the current check if it's high enough to
1438  // trigger a random change in ROn/ROff
1439  // and this is the start of a new step, then get a new
1440  // gaussian distributed variable
1441  if( ri.model_.randomResNoiseOn_ &&
1444  {
1447  // get a poisson distributed number with mean lambda
1448  //int aNum = ri.model_.randomNumberGen_->poissonRandom(ri.model_.randomResNoiseLambda_);
1449  // the actual update interval is the poisson distributed number times the base update time.
1450  //ri.resNoiseNextUpdateTime = aNum * ri.model_.randomResUpdateTime_;
1451  // don't allow a zero time for the next update.
1452  //if( aNum == 0 )
1453  // ri.resNoiseNextUpdateTime = ri.model_.randomResEpsilonUpdateTime_;
1455 
1456  // update the resistance factor
1457  if( ri.resNoiseHiLoState == 0 )
1458  {
1459  // in low state so move to high state
1460  ri.resNoiseHiLoState = 1;
1462  }
1463  else
1464  {
1465  // in high state so move to low state
1466  ri.resNoiseHiLoState = 0;
1468  }
1469 
1470  // ri.resNoiseRfactor = ri.model_.randomNumberGen_->gaussianRandom(ri.model_.randomResNoiseMean_, ri.model_.randomResNoiseSD_);
1471  // Xyce::dout() << " updating random component of R..." << ri.resNoiseRfactor << " Reff was " << ri.Reff;
1472  //Xyce::dout() << " updating random component of R: p( " << ri.model_.randomResNoiseLambda_<< " ) = " << 0
1473  // << " and next update time interval is " << ri.resNoiseNextUpdateTime << std::endl;
1474  {
1475  // extra scope to limit variable extent. Otherwise one will have lots of
1476  // similar variables
1477  Sacado::Fad::SFad<double,1> varX( 1, 0, x );
1478  Sacado::Fad::SFad<double,1> paramRon( ri.resNoiseRfactor*ri.model_.ROn_ );
1479  Sacado::Fad::SFad<double,1> paramRoff(ri.resNoiseRfactor*ri.model_.ROff_ );
1480  Sacado::Fad::SFad<double,1> paramXon( ri.model_.xOn_ * ri.model_.xScaling_ );
1481  Sacado::Fad::SFad<double,1> paramXoff( ri.model_.xOff_ * ri.model_.xScaling_ );
1482  Sacado::Fad::SFad<double,1> resultFad;
1483  if(ri.IVrelation==0)
1484  {
1485  // liner current voltage relationship
1486  //Xyce::dout() << ri.getName() << "About to call ReffLin..." << ri.model_.ROn_ << " , " << ri.model_.ROff_ << ", " << ri.model_.xOn_ << ", " << ri.model_.xOff_ << std::endl;
1487 
1488  ReffLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
1489  //fval = Ron + (X - Xon) * (Roff-Ron)/(Xoff - Xon);
1490  //dReff/dV1 = 0
1491  //dReff/dV2 = 0;
1492  //dReff/dx = (Roff-Ron)/(Xoff - Xon);
1493  }
1494  else if( ri.IVrelation==1 )
1495  {
1496  //Xyce::dout() << ri.getName() << "About to call ReffNonLin..." << std::endl;
1497  ReffNonLin( varX, paramRon, paramRoff, paramXon, paramXoff, resultFad );
1498  }
1499  ri.Reff = resultFad.val();
1500  ri.dReffdx = resultFad.dx(0);
1501  ri.dReffdvpos = (v_pos-v_neg)*(-pow( ri.Reff, 2 ) * ri.dReffdx) / (ri.model_.xScaling_);
1502  ri.dReffdvneg = -(v_pos-v_neg)*(-pow( ri.Reff, 2 ) * ri.dReffdx) / (ri.model_.xScaling_);
1503  //ri.dReffdvpos = 0.0;
1504  //ri.dReffdvneg = 0.0;
1505 
1506  }
1507  //Xyce::dout() << " Reff now " << ri.Reff << std::endl;
1508  ri.G = 1.0/(ri.Reff);
1509  ri.i0 = (v_pos-v_neg)*ri.G;
1510  }
1511 
1512  if( ri.model_.randomResNoiseOn_ )
1513  {
1514  stoVec[ri.li_store_tdt] = ri.resNoiseNextUpdateTime;
1515  }
1516 
1517  {
1518  // calculate the xVar F equation update and the
1519  // derivatives d(F(x))/dVpos, d(F(x))/dVneg and d(F(x))/dx
1520  Sacado::Fad::SFad<double,3> varVpos( 3, 0, v_pos);
1521  Sacado::Fad::SFad<double,3> varVneg( 3, 1, v_neg);
1522  Sacado::Fad::SFad<double,3> varX( 3, 2, x );
1523  Sacado::Fad::SFad<double,3> parami( ri.i0 );
1524  Sacado::Fad::SFad<double,3> paramG( ri.G );
1525  Sacado::Fad::SFad<double,3> paramiOff( ri.model_.iOff_ );
1526  Sacado::Fad::SFad<double,3> paramiOn( ri.model_.iOn_ );
1527 
1528  //Sacado::Fad::SFad<double,3> paramkOff( ri.model_.kOff_ );
1529  //Sacado::Fad::SFad<double,3> paramkOn( ri.model_.kOn_ );
1530  Sacado::Fad::SFad<double,3> paramkOff( ri.model_.kOff_*ri.model_.xScaling_ );
1531  Sacado::Fad::SFad<double,3> paramkOn( ri.model_.kOn_ *ri.model_.xScaling_ );
1532  Sacado::Fad::SFad<double,3> paramAlphaOff( ri.model_.alphaOff_ );
1533  Sacado::Fad::SFad<double,3> paramAlphaOn( ri.model_.alphaOn_ );
1534  Sacado::Fad::SFad<double,3> resultFad;
1535 
1536  xVarFterm(varVpos, varVneg, varX, paramG, paramiOff, paramiOn, paramkOff, paramkOn,
1537  paramAlphaOff, paramAlphaOn, resultFad);
1538 
1539  // calculate the window function value
1540  Sacado::Fad::SFad<double,3> windowFunctionFad;
1541 
1542  if( ri.model_.windowType_ == 0)
1543  {
1544  // no window
1545  windowFunctionFad = 1.0;
1546  }
1547  else if( ri.model_.windowType_ == 1)
1548  {
1549  Sacado::Fad::SFad<double,3> paramD( ri.model_.D_ );
1550  Sacado::Fad::SFad<double,3> paramP( ri.model_.p_ );
1551  JogelkarWindowFunction( varX, paramD, paramP, windowFunctionFad );
1552  }
1553  else if( ri.model_.windowType_ == 2)
1554  {
1555  Sacado::Fad::SFad<double,3> paramD( ri.model_.D_ );
1556  Sacado::Fad::SFad<double,3> paramP( ri.model_.p_ );
1557  BiolekWindowFunction( varX, paramD, paramP, parami, windowFunctionFad );
1558  }
1559  else if( ri.model_.windowType_ == 3)
1560  {
1561  Sacado::Fad::SFad<double,3> paramD( ri.model_.D_ );
1562  Sacado::Fad::SFad<double,3> paramP( ri.model_.p_ );
1563  Sacado::Fad::SFad<double,3> paramJ( ri.model_.j_ );
1564  ProdromakisWindowFunction( varX, paramD, paramP, paramJ, windowFunctionFad );
1565  }
1566  else if( ri.model_.windowType_ == 4)
1567  {
1568  Sacado::Fad::SFad<double,3> paramAOn( ri.model_.aOn_ *ri.model_.xScaling_ );
1569  Sacado::Fad::SFad<double,3> paramAOff( ri.model_.aOff_ *ri.model_.xScaling_ );
1570  Sacado::Fad::SFad<double,3> paramWc( ri.model_.wc_ *ri.model_.xScaling_ );
1571  TEAMWindowFunctionF( varX, parami, paramAOff, paramAOn, paramWc, windowFunctionFad );
1572  }
1573 
1574  // xVarFterm * WindowFunction
1575  // dertivative terms then become: WindowFunction dxVarFterm/dv1 + xVarFterm dW/dv1
1576 
1577  ri.xVarFContribution = resultFad.val()*windowFunctionFad.val();
1578  ri.dxFEqdVpos = windowFunctionFad.val()*resultFad.dx(0);
1579  ri.dxFEqdVneg = windowFunctionFad.val()*resultFad.dx(1);
1580  ri.dxFEqdx = windowFunctionFad.val()*resultFad.dx(2) + resultFad.val()*windowFunctionFad.dx(2) / (ri.model_.xScaling_);
1581 
1582  }
1583 
1584  //ri.G = 1.0/(rfactor*ri.Reff);
1585  //ri.i0 = (v_pos-v_neg)*ri.G;
1586 
1587  }
1588 
1589  return true;
1590 }
1591 
1592 //-----------------------------------------------------------------------------
1593 // Function : Xyce::Device::MemristorTEAM::Master::loadDAEVectors
1594 // Purpose :
1595 // Special Notes :
1596 // Scope : public
1597 // Creator : Richard Schiek, Electrical Models & Simulations
1598 // Creation Date : 10/23/2014
1599 //-----------------------------------------------------------------------------
1600 //
1601 // Load DAE vectors of all memristor instances, regardless of model
1602 //
1603 // @param solVec solution vector
1604 // @param fVec f vector
1605 // @param qVec q vector
1606 // @param storeLeadF store lead current f vector
1607 // @param storeLeadQ store lead current q vector
1608 //
1609 // @return true on success
1610 //
1611 // @note Because the memristor device re-implements the base-class
1612 // Master::loadDAEVectors, the Instance::loadDAEFVector method is
1613 // never called. This method replaces those, and does the same work
1614 // but inside a loop over all memristor instances.
1615 //
1616 // @see Xyce::Device::MemristorTEAM::Instance::loadDAEFVector
1617 //
1618 // @author Eric Keiter, SNL
1619 // @date 11/26/08
1620 
1621 bool Master::loadDAEVectors (double * solVec, double * fVec, double *qVec, double * bVec, double * storeLeadF, double * storeLeadQ, double * leadF, double * leadQ, double * junctionV)
1622 {
1623  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1624  {
1625  Instance & ri = *(*it);
1626  fVec[ri.li_Pos] += ri.i0;
1627  fVec[ri.li_Neg] += -ri.i0;
1628  fVec[ri.li_x] += ri.xVarFContribution;
1629  qVec[ri.li_x] += solVec[ri.li_x];
1630  storeLeadF[ri.li_store_R] = 1.0/ri.G;
1631 
1632  if( ri.loadLeadCurrent )
1633  {
1634  leadF[ri.li_branch_data] = ri.i0;
1635  junctionV[ri.li_branch_data] = solVec[ri.li_Pos] - solVec[ri.li_Neg];
1636  }
1637 
1638  }
1639  return true;
1640 }
1641 
1642 //-----------------------------------------------------------------------------
1643 // Function : Xyce::Device::MemristorTEAM::Master::loadDAEMatrices
1644 // Purpose :
1645 // Special Notes :
1646 // Scope : public
1647 // Creator : Richard Schiek, Electrical Models & Simulations
1648 // Creation Date : 10/23/2014
1649 //-----------------------------------------------------------------------------
1650 //
1651 // Load DAE matrices for all memristor instances, regardless of model
1652 //
1653 // @param dFdx matrix of derivatives of F vector with respect to solution
1654 // @param dQdx matrix of derivatives of Q vector with respect to solution
1655 //
1656 // @return true on success
1657 //
1658 // @note Because the memristor device re-implements the base-class
1659 // Master::loadDAEMatrices, the Instance::loadDAEdFdx method is
1660 // never called. This method replaces those, and does the same work
1661 // but inside a loop over all memristor instances.
1662 //
1663 // @see Xyce::Device::MemristorTEAM::Instance::loadDAEdFdx
1664 //
1665 // @author Eric Keiter, SNL
1666 // @date 11/26/08
1667 
1668 bool Master::loadDAEMatrices(Linear::Matrix & dFdx, Linear::Matrix & dQdx)
1669 {
1670  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1671  {
1672  Instance & ri = *(*it);
1673 
1674 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
1675  *(ri.f_PosEquPosNodePtr) += ri.G;
1676  *(ri.f_PosEquNegNodePtr) -= ri.G;
1677  *(ri.f_PosEquXNodePtr) += ri.dReffdvpos;
1678  *(ri.f_NegEquPosNodePtr) -= ri.G;
1679  *(ri.f_NegEquNegNodePtr) += ri.G;
1680  *(ri.f_NegEquXNodePtr ) += ri.dReffdvneg;
1681  *(ri.f_XEquPosNodePtr ) += ri.dxFEqdVpos;
1682  *(ri.f_XEquNegNodePtr ) += ri.dxFEqdVneg;
1683  *(ri.f_XEquXNodePtr ) += ri.dxFEqdx;
1684  *(ri.q_XEquXNodePtr ) = 1.0;
1685 #else
1686  dFdx[ri.li_Pos][ri.APosEquPosNodeOffset] += ri.G;
1687  dFdx[ri.li_Pos][ri.APosEquNegNodeOffset] -= ri.G;
1688  dFdx[ri.li_Pos][ri.APosEquXNodeOffset] += ri.dReffdvpos;
1689  dFdx[ri.li_Neg][ri.ANegEquPosNodeOffset] -= ri.G;
1690  dFdx[ri.li_Neg][ri.ANegEquNegNodeOffset] += ri.G;
1691  dFdx[ri.li_Neg][ri.ANegEquXNodeOffset] += ri.dReffdvneg;
1692  dFdx[ri.li_x][ri.XEquVPosOffset] += ri.dxFEqdVpos;
1693  dFdx[ri.li_x][ri.XEquVNegOffset] += ri.dxFEqdVneg;
1694  dFdx[ri.li_x][ri.XEquXOffset] += ri.dxFEqdx;
1695  dQdx[ri.li_x][ri.XEquXOffset] = 1.0;
1696 #endif
1697  }
1698 
1699  return true;
1700 }
1701 
1702 //-----------------------------------------------------------------------------
1703 // Function : Xyce::Device::MemristorTEAM::Traits::factory
1704 // Purpose :
1705 // Special Notes :
1706 // Scope : public
1707 // Creator : Richard Schiek, Electrical Models & Simulations
1708 // Creation Date : 10/23/2014
1709 //-----------------------------------------------------------------------------
1710 //
1711 // Create a new instance of the MemristorTEAM device.
1712 //
1713 // @param configuration
1714 // @param factory_block
1715 //
1716 Device *
1717 Traits::factory(const Configuration &configuration, const FactoryBlock &factory_block)
1718 {
1719  return new Master(configuration, factory_block, factory_block.solverState_, factory_block.deviceOptions_);
1720 }
1721 
1722 //-----------------------------------------------------------------------------
1723 // Function : Xyce::Device::MemristorTEAM::registerDevice
1724 // Purpose :
1725 // Special Notes :
1726 // Scope : public
1727 // Creator : Richard Schiek, Electrical Models & Simulations
1728 // Creation Date : 10/23/2014
1729 //-----------------------------------------------------------------------------
1730 //
1731 // Define how to use the device in a netlist.
1732 //
1733 // This method is called from the Xyce::Device::registerOpenDevices
1734 // function, which in turn is called by the device manager.
1735 //
1736 // The device is declared here to be an "memristor" device, which must
1737 // take a model card of type "memristor". This device will correspond to model
1738 // level 2 of memristor models.
1739 void
1741 {
1743  .registerDevice("memristor", 2)
1744  .registerModelType("memristor", 2);
1745 }
1746 
1747 //-----------------------------------------------------------------------------
1748 // Function : memristorTEAMSensitivity::operator
1749 // Purpose : produces df/dp and dq/dp, where p=R.
1750 // Special Notes :
1751 // Scope : public
1752 // Creator : Richard Schiek, Electrical Models & Simulations
1753 // Creation Date : 10/23/2014
1754 //-----------------------------------------------------------------------------
1756  const ParameterBase &entity,
1757  const std::string & name,
1758  std::vector<double> & dfdp,
1759  std::vector<double> & dqdp,
1760  std::vector<double> & dbdp,
1761  std::vector<int> & Findices,
1762  std::vector<int> & Qindices,
1763  std::vector<int> & Bindices
1764  ) const
1765 {
1766  const ParameterBase * e1 = &entity;
1767  const Instance * in = dynamic_cast<const Instance *> (e1);
1768 
1769  double * solVec = in->extData.nextSolVectorRawPtr;
1770  double v_pos = solVec[in->li_Pos];
1771  double v_neg = solVec[in->li_Neg];
1772 
1773  double dfdpLoc = -(v_pos-v_neg)*in->G*in->G;
1774 
1775  dfdp.resize(2);
1776  dfdp[0] = +dfdpLoc;
1777  dfdp[1] = -dfdpLoc;
1778 
1779  Findices.resize(2);
1780  Findices[0] = in->li_Pos;
1781  Findices[1] = in->li_Neg;
1782 }
1783 
1784 } // namespace MemristorTEAM
1785 } // namespace Device
1786 } // namespace Xyce
const InstanceName & getName() const
const DeviceOptions & deviceOptions_
Descriptor & addPar(const char *parName, T default_value, T U::*varPtr)
Adds the parameter description to the parameter map.
Definition: N_DEV_Pars.h:1429
virtual void registerStoreLIDs(const std::vector< int > &stoLIDVecRef)
virtual bool loadDAEMatrices(Linear::Matrix &dFdx, Linear::Matrix &dQdx)
Populates the device's Jacobian object with these pointers.
void G(const ScalarT &V1, const ScalarT &V2, const ScalarT &Ap, const ScalarT &An, ScalarT &Vp, ScalarT &Vn, ScalarT &fval)
Pure virtual class to augment a linear system.
void addInternalNode(Util::SymbolTable &symbol_table, int index, const InstanceName &instance_name, const std::string &lead_name)
virtual void registerLIDs(const std::vector< int > &intLIDVecRef, const std::vector< int > &extLIDVecRef)
void setNumStoreVars(int num_store_vars)
void addBranchDataNode(Util::SymbolTable &symbol_table, int index, const InstanceName &instance_name, const std::string &lead_name)
InstanceVector::const_iterator getInstanceEnd() const
Returns an iterator to the ending of the vector of all instances created for this device...
Base class for all parameters.
Definition: N_DEV_Pars.h:169
virtual std::ostream & printOutInstances(std::ostream &os) const
#define AssertLIDs(cmp)
virtual void registerJacLIDs(const JacobianStamp &jacLIDVec)
static void loadInstanceParameters(ParametricData< Instance > &p)
virtual void operator()(const ParameterBase &entity, const std::string &name, std::vector< double > &dfdp, std::vector< double > &dqdp, std::vector< double > &dbdp, std::vector< int > &Findices, std::vector< int > &Qindices, std::vector< int > &Bindices) const
InstanceVector::const_iterator getInstanceBegin() const
Returns an iterator to the beginning of the vector of all instances created for this device...
std::vector< Param > params
Parameters from the line.
void JogelkarWindowFunction(const ScalarT &w, const ScalarT &D, const ScalarT &p, ScalarT &fval)
void setParams(const std::vector< Param > &params)
const std::string & getName() const
The FactoryBlock contains parameters needed by the device, instance and model creation functions...
virtual void registerJacLIDs(const std::vector< std::vector< int > > &jacLIDVec)
void xVarFterm(const ScalarT &Vpos, const ScalarT &Vneg, const ScalarT &x, const ScalarT &G, const ScalarT &iOff, const ScalarT &iOn, const ScalarT &kOff, const ScalarT &kOn, const ScalarT &alphaOff, const ScalarT &alphaOn, ScalarT &fval)
Instance(const Configuration &configuration, const InstanceBlock &instance_block, Model &model, const FactoryBlock &factory_block)
void ProdromakisWindowFunction(const ScalarT &w, const ScalarT &D, const ScalarT &p, const ScalarT &j, ScalarT &fval)
static Config< T > & addConfiguration()
Adds the device to the Xyce device configuration.
Linear::Matrix * dFdxMatrixPtr
void ReffLin(const ScalarT &X, const ScalarT &Ron, const ScalarT &Roff, const ScalarT &Xon, const ScalarT &Xoff, ScalarT &fval)
virtual bool processInstanceParams()
processInstanceParams
The Device class is an interface for device implementations.
Definition: N_DEV_Device.h:101
void addStoreNode(Util::SymbolTable &symbol_table, int index, const InstanceName &instance_name, const std::string &lead_name)
virtual bool updateTemperature(const double &temp_tmp)
int timeStepNumber_
Memristor, LTRA, TRA, testing if debug or jacobian for testing.
const SolverState & solverState_
Class Configuration contains device configuration data.
virtual bool loadDAEVectors(double *solVec, double *fVec, double *qVec, double *bVec, double *storeLeadF, double *storeLeadQ, double *leadF, double *leadQ, double *junctionV)
Populates the device's ExternData object with these pointers.
virtual void registerBranchDataLIDs(const std::vector< int > &branchLIDVecRef)
virtual bool updateState(double *solVec, double *staVec, double *stoVec)
Updates the devices state information.
virtual void forEachInstance(DeviceInstanceOp &op) const
const SolverState & getSolverState() const
virtual void registerStateLIDs(const std::vector< int > &staLIDVecRef)
void setNumBranchDataVars(int num_branch_data_vars)
void ReffNonLin(const ScalarT &X, const ScalarT &Ron, const ScalarT &Roff, const ScalarT &Xon, const ScalarT &Xoff, ScalarT &fval)
static void loadModelParameters(ParametricData< Model > &p)
#define Xyce_NONPOINTER_MATRIX_LOAD
Definition: N_DEV_Bsrc.C:97
static std::vector< std::vector< int > > jacStamp
void loadNodeSymbols(Util::SymbolTable &symbol_table) const
Populates and returns the store name map.
ModelBlock represents a .MODEL line from the netlist.
Manages parameter binding for class C.
Definition: N_DEV_Pars.h:214
InstanceBlock represent a device instance line from the netlist.
double currTime_
DeviceEntity for expression time, breakpoints DeviceMgr for dependent parameters, breakpoints...
std::vector< Param > params
Xyce::Util::RandomNumbers * randomNumberGen_
Linear::Matrix * dQdxMatrixPtr
static Device * factory(const Configuration &configuration, const FactoryBlock &factory_block)
void BiolekWindowFunction(const ScalarT &w, const ScalarT &D, const ScalarT &p, const ScalarT &i, ScalarT &fval)
void TEAMWindowFunctionF(const ScalarT &x, const ScalarT &i, const ScalarT &aOff, ScalarT &aOn, const ScalarT &wc, ScalarT &fval)
virtual bool processParams()
processParams
const SolverState & getSolverState() const
Returns the solver state given during device construction.
void setModParams(const std::vector< Param > &params)