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