Xyce  6.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
N_DEV_Inductor.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------------
2 // Copyright Notice
3 //
4 // Copyright 2002 Sandia Corporation. Under the terms
5 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S.
6 // Government retains certain rights in this software.
7 //
8 // Xyce(TM) Parallel Electrical Simulator
9 // Copyright (C) 2002-2014 Sandia Corporation
10 //
11 // This program is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // This program is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU General Public License
22 // along with this program. If not, see <http://www.gnu.org/licenses/>.
23 //-----------------------------------------------------------------------------
24 
25 //-------------------------------------------------------------------------
26 // Filename : $RCSfile: N_DEV_Inductor.C,v $
27 //
28 // Purpose :
29 //
30 // Special Notes :
31 //
32 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
33 //
34 // Creation Date : 02/28/00
35 //
36 // Revision Information:
37 // ---------------------
38 //
39 // Revision Number: $Revision: 1.273.2.1 $
40 //
41 // Revision Date : $Date: 2014/08/28 18:24:58 $
42 //
43 // Current Owner : $Author: tvrusso $
44 //-------------------------------------------------------------------------
45 #include <Xyce_config.h>
46 
47 // ---------- Standard Includes ----------
48 #include <N_UTL_Misc.h>
49 
50 #ifdef HAVE_CMATH
51 #include <cmath>
52 #else
53 #include <math.h>
54 #endif
55 
56 #include <algorithm>
57 #include <fstream>
58 #include <set>
59 #include <vector>
60 
61 // ---------- Xyce Includes ----------
62 #include <N_DEV_DeviceOptions.h>
63 #include <N_DEV_ExternData.h>
64 #include <N_DEV_Inductor.h>
65 #include <N_DEV_MatrixLoadData.h>
66 #include <N_DEV_SolverState.h>
67 #include <N_DEV_Message.h>
68 #include <N_ERH_ErrorMgr.h>
69 
70 #include <N_LAS_Vector.h>
71 #include <N_LAS_Matrix.h>
72 
73 namespace Xyce {
74 namespace Device {
75 namespace Inductor {
76 
78 {
79  p.addPar("L", 0.0, &Inductor::Instance::baseL)
80  .setExpressionAccess(ParameterType::TIME_DEP)
81  .setUnit(U_HENRY)
82  .setDescription("Inductance")
83  .setAnalyticSensitivityAvailable(true)
84  .setSensitivityFunctor(&indSens);
85 
86  p.addPar("IC", 0.0, &Inductor::Instance::IC)
87  .setGivenMember(&Inductor::Instance::ICGiven)
88  .setUnit(U_AMP)
89  .setDescription("Initial current through device");
90 
91  p.addPar("TEMP", 0.0, &Inductor::Instance::temp)
92  .setExpressionAccess(ParameterType::TIME_DEP)
93  .setGivenMember(&Inductor::Instance::tempGiven)
94  .setUnit(U_DEGC)
95  .setCategory(CAT_MATERIAL)
96  .setDescription("Device temperature");
97 
99  .setGivenMember(&Inductor::Instance::tempCoeff1Given)
100  .setUnit(U_DEGCM1)
101  .setDescription("Linear Temperature Coefficient");
102 
103  p.addPar("TC2", 0.0, &Inductor::Instance::tempCoeff2)
104  .setGivenMember(&Inductor::Instance::tempCoeff2Given)
105  .setUnit(U_DEGCM2)
106  .setDescription("Quadratic Temperature Coefficient");
107 
108  // This call tells the parameter handling code that TC can be specified
109  // as a vector with up to two elements as in TC=a,b. It then translates
110  // TC=a,b into TC1=a TC2=b. Likewise, TC=a will translate into TC1=a
111  p.makeVector ("TC", 2);
112 }
113 
115 {
116  // Set up double precision variables:
118  .setUnit(U_NONE)
119  .setDescription("Inductance Multiplier");
120 
121  p.addPar("IC", 0.0, &Inductor::Model::IC)
122  .setUnit(U_AMP)
123  .setDescription("Initial current through device");
124 
125  p.addPar("TNOM", 27.0, &Inductor::Model::tnom)
126  .setUnit(U_DEGC)
127  .setCategory(CAT_MATERIAL)
128  .setDescription("Reference temperature");
129 
130  p.addPar("TC1",0.0, &Inductor::Model::tempCoeff1)
131  .setUnit(U_DEGCM1)
132  .setCategory(CAT_MATERIAL)
133  .setDescription("First order temperature coeff.");
134 
135  p.addPar("TC2", 0.0, &Inductor::Model::tempCoeff2)
136  .setUnit(U_DEGCM2)
137  .setCategory(CAT_MATERIAL)
138  .setDescription("Second order temperature coeff.");
139 }
140 
141 //
142 // static class member inits
143 //
144 std::vector< std::vector<int> > Instance::jacStamp_BASE;
145 
146 // Class Instance
147 //-----------------------------------------------------------------------------
148 // Function : Instance::processParams
149 // Purpose :
150 // Special Notes :
151 // Scope : public
152 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
153 // Creation Date : 6/03/02
154 //-----------------------------------------------------------------------------
156 {
157  // If there are any time dependent parameters, set their values at for
158  // the current time.
159 
160  // now set the temperature related stuff.
162  return true;
163 }
164 
165 //-----------------------------------------------------------------------------
166 // Function : Instance::updateTemperature
167 // Purpose :
168 // Special Notes :
169 // Scope : public
170 // Creator : Tom Russo, Component Information and Models
171 // Creation Date : 02/27/01
172 //-----------------------------------------------------------------------------
173 bool Instance::updateTemperature ( const double & temp)
174 {
175  double difference = temp - model_.tnom;
176  double factor = model_.inductanceMultiplier*(1.0 + tempCoeff1*difference +
177  tempCoeff2*difference*difference);
178  L = baseL*factor;
179  return true;
180 }
181 
182 //-----------------------------------------------------------------------------
183 // Function : Instance::Instance
184 // Purpose : constructor
185 // Special Notes :
186 // Scope : public
187 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
188 // Creation Date : 3/16/00
189 //-----------------------------------------------------------------------------
191  const Configuration & configuration,
192  const InstanceBlock & instance_block,
193  Model & model,
194  const FactoryBlock & factory_block)
195  : DeviceInstance(instance_block, configuration.getInstanceParameters(), factory_block),
196  L(0),
197  IC(0),
198  ICGiven(false),
199  model_(model),
200  baseL(0.0),
201  temp(getDeviceOptions().temp.getImmutableValue<double>()),
202  tempGiven(0),
203  tempCoeff1(0.0),
204  tempCoeff2(0.0),
205  tempCoeff1Given(false),
206  tempCoeff2Given(false),
207  li_fstate(-1),
208  li_Pos(-1),
209  li_Neg(-1),
210  li_Bra(-1),
211  ABraEquPosNodeOffset(-1),
212  ABraEquNegNodeOffset(-1),
213  ABraEquBraVarOffset(-1),
214  APosEquBraVarOffset(-1),
215  ANegEquBraVarOffset(-1)
217  ,
218  fPosEquBraVarPtr(0),
219  fNegEquBraVarPtr(0),
220  fBraEquPosNodePtr(0),
221  fBraEquNegNodePtr(0),
222  fBraEquBraVarPtr(0),
223  qBraEquBraVarPtr(0)
224 #endif
225 {
226  numExtVars = 2;
227  numIntVars = 1;
228  numStateVars = 1;
229 
230  if( jacStamp_BASE.empty() )
231  {
232  jacStamp_BASE.resize(3);
233 
234  jacStamp_BASE[0].resize(1);
235  jacStamp_BASE[0][0] = 2;
236 
237  jacStamp_BASE[1].resize(1);
238  jacStamp_BASE[1][0] = 2;
239 
240  jacStamp_BASE[2].resize(3);
241  jacStamp_BASE[2][0] = 0;
242  jacStamp_BASE[2][1] = 1;
243  jacStamp_BASE[2][2] = 2;
244  }
245 
246  // Set params to constant default values:
247  setDefaultParams ();
248 
249  // Set params according to instance line and constant defaults from metadata:
250  setParams (instance_block.params);
251 
252  // Set any non-constant parameter defaults:
253  if (!given("L"))
254  {
255  UserError0(*this) << "Could not find L parameter in instance.";
256  }
257 
258  if (!given("TEMP"))
259  temp = getDeviceOptions().temp.getImmutableValue<double>();
260 
261  if (!tempCoeff1Given)
263  if (!tempCoeff2Given)
265 
266  // Calculate any parameters specified as expressions:
268 
269  // calculate dependent (ie computed) params and check for errors:
270  processParams ();
271 
272  // set up numIntVars:
273  numIntVars = 1;
274 
275  // set up numStateVars:
276  numStateVars = 2;
277 }
278 
279 //-----------------------------------------------------------------------------
280 // Function : Instance::~Instance
281 // Purpose : destructor
282 // Special Notes :
283 // Scope : public
284 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
285 // Creation Date : 3/16/00
286 //-----------------------------------------------------------------------------
288 {
289 }
290 
291 //-----------------------------------------------------------------------------
292 // Function : Instance::setupPointers
293 // Purpose :
294 // Special Notes :
295 // Scope : public
296 // Creator : Eric Keiter, SNL
297 // Creation Date : 11/30/08
298 //-----------------------------------------------------------------------------
300 {
301 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
302  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
303  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
304 
310 
312 
313 #endif
314 }
315 
316 //-----------------------------------------------------------------------------
317 // Function : Instance::registerLIDs
318 // Purpose :
319 // Special Notes :
320 // Scope : public
321 // Creator : Robert Hoekstra, SNL, Parallel Computational Sciences
322 // Creation Date : 6/21/02
323 //-----------------------------------------------------------------------------
324 void Instance::registerLIDs(const std::vector<int> & intLIDVecRef,
325  const std::vector<int> & extLIDVecRef)
326 {
327  AssertLIDs(intLIDVecRef.size() == numIntVars);
328  AssertLIDs(extLIDVecRef.size() == numExtVars);
329 
330  // copy over the global ID lists.
331  intLIDVec = intLIDVecRef;
332  extLIDVec = extLIDVecRef;
333 
334  // Now use these lists to obtain the indices into the
335  // linear algebra entities. This assumes an order.
336  // For the matrix indices, first do the rows.
337 
338  li_Pos = extLIDVec[0];
339  li_Neg = extLIDVec[1];
340  li_Bra = intLIDVec[0];
341 
342 #ifdef Xyce_DEBUG_DEVICE
343  if (getDeviceOptions().debugLevel > 0)
344  {
345  Xyce::dout() << section_divider << std::endl;
346 
347  Xyce::dout() << "::registerLIDs:\n";
348  Xyce::dout() << " name = " << getName() << std::endl;
349 
350  Xyce::dout() << "\nlocal solution indices:\n";
351  Xyce::dout() << " li_Pos = "<< li_Pos << std::endl;
352  Xyce::dout() << " li_Neg = "<< li_Neg << std::endl;
353  Xyce::dout() << " li_Bra = "<< li_Bra << std::endl;
354 
355  Xyce::dout() << section_divider << std::endl;
356  }
357 #endif
358 
359 }
360 
361 //-----------------------------------------------------------------------------
362 // Function : Instance::getIntNameMap
363 // Purpose :
364 // Special Notes :
365 // Scope : public
366 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
367 // Creation Date : 05/13/05
368 //-----------------------------------------------------------------------------
369 std::map<int,std::string> & Instance::getIntNameMap ()
370 {
371  // set up the internal name map, if it hasn't been already.
372  if (intNameMap.empty ())
373  {
374  intNameMap[ li_Bra ] = spiceInternalName(getName(), "branch");
375  }
376 
377  return intNameMap;
378 }
379 
380 //-----------------------------------------------------------------------------
381 // Function : Instance::registerStateLIDs
382 // Purpose :
383 // Special Notes :
384 // Scope : public
385 // Creator : Robert Hoekstra, SNL, Parallel Computational Sciences
386 // Creation Date : 6/22/02
387 //-----------------------------------------------------------------------------
388 void Instance::registerStateLIDs( const std::vector<int> & staLIDVecRef )
389 {
390  AssertLIDs(staLIDVecRef.size() == numStateVars);
391 
392  // copy over the global ID lists.
393  staLIDVec = staLIDVecRef;
394 
395  li_fstate = staLIDVec[0];
396 }
397 
398 //-----------------------------------------------------------------------------
399 // Function : Instance::jacobianStamp
400 // Purpose :
401 // Special Notes :
402 // Scope : public
403 // Creator : Robert Hoekstra, SNL
404 // Creation Date : 08/21/02
405 //-----------------------------------------------------------------------------
406 const std::vector< std::vector<int> > & Instance::jacobianStamp() const
407 {
408  return jacStamp_BASE;
409 }
410 
411 //-----------------------------------------------------------------------------
412 // Function : Instance::registerJacLIDs
413 // Purpose :
414 // Special Notes :
415 // Scope : public
416 // Creator : Robert Hoekstra, SNL
417 // Creation Date : 08/27/02
418 //-----------------------------------------------------------------------------
419 void Instance::registerJacLIDs( const std::vector< std::vector<int> > & jacLIDVec )
420 {
421  DeviceInstance::registerJacLIDs( jacLIDVec );
422 
423  APosEquBraVarOffset = jacLIDVec[0][0];
424  ANegEquBraVarOffset = jacLIDVec[1][0];
425  ABraEquPosNodeOffset = jacLIDVec[2][0];
426  ABraEquNegNodeOffset = jacLIDVec[2][1];
427  ABraEquBraVarOffset = jacLIDVec[2][2];
428 }
429 
430 //-----------------------------------------------------------------------------
431 // Function : Instance::updatePrimaryState
432 // Purpose :
433 // Special Notes :
434 // Scope : public
435 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
436 // Creation Date : 02/01/01
437 //-----------------------------------------------------------------------------
439 {
440  double * solVec = extData.nextSolVectorRawPtr;
441  double * staVec = extData.nextStaVectorRawPtr;
442  double current = solVec[li_Bra];
443  if( (getSolverState().dcopFlag) && ICGiven )
444  current = IC;
445 
446  f0 = L*current;
447  staVec[li_fstate] = f0;
448 
449  return true;
450 }
451 
452 //-----------------------------------------------------------------------------
453 // Function : Instance::updateSecondaryState
454 // Purpose :
455 // Special Notes :
456 // Scope : public
457 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
458 // Creation Date : 02/01/01
459 //-----------------------------------------------------------------------------
461 {
462  return true;
463 }
464 
465 //-----------------------------------------------------------------------------
466 // Function : Instance::loadDAEQVector
467 //
468 // Purpose : Loads the Q-vector contributions for a single
469 // instance.
470 //
471 // Special Notes : The "Q" vector is part of a standard DAE formalism in
472 // which the system of equations is represented as:
473 //
474 // f(x) = dQ(x)/dt + F(x) - B(t) = 0
475 //
476 // Scope : public
477 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
478 // Creation Date : 01/26/03
479 //-----------------------------------------------------------------------------
481 {
482  double * qVec = extData.daeQVectorRawPtr;
483 
484  qVec[li_Bra] += f0;
485 
486  return true;
487 }
488 
489 //-----------------------------------------------------------------------------
490 // Function : Instance::loadDAEFVector
491 //
492 // Purpose : Loads the F-vector contributions for a single
493 // instance.
494 //
495 // Special Notes :
496 //
497 // Scope : public
498 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
499 // Creation Date : 01/26/03
500 //-----------------------------------------------------------------------------
502 {
503  double current;
504  double coef;
505 
506  double * solVec = extData.nextSolVectorRawPtr;
507  double * fVec = extData.daeFVectorRawPtr;
508 
509  double v_pos = solVec[li_Pos];
510  double v_neg = solVec[li_Neg];
511  double vind = v_pos-v_neg;
512 
513  // In the case that an initial condition is specified, the inductor is set up
514  // like a current source for the DC operating point. We don't deal with the
515  // node voltages in that case, so set coef to 0.
516  if (getSolverState().dcopFlag && ICGiven)
517  {
518  current = IC;
519  coef = 0.0;
520  }
521  else
522  {
523  current = solVec[li_Bra];
524  coef = -vind;
525  }
526 
527  // load the current into the two KCL rhs vector elements
528  fVec[li_Pos] += current;
529  fVec[li_Neg] += -current;
530  fVec[li_Bra] += coef;
531 
532  return true;
533 }
534 
535 //-----------------------------------------------------------------------------
536 // Function : Instance::loadDAEdQdx
537 //
538 // Purpose : Loads the Q-vector contributions for a single
539 // instance.
540 // Scope : public
541 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
542 // Creation Date : 03/05/04
543 //-----------------------------------------------------------------------------
545 {
546  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
547  dQdx[li_Bra][ABraEquBraVarOffset] += L;
548  return true;
549 }
550 
551 //-----------------------------------------------------------------------------
552 // Function : Instance::loadDAEdFdx ()
553 //
554 // Purpose : Loads the F-vector contributions for a single
555 // instance.
556 //
557 // Special Notes :
558 //
559 // Scope : public
560 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
561 // Creation Date : 03/05/04
562 //-----------------------------------------------------------------------------
564 {
565  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
566  if ( getSolverState().dcopFlag && ICGiven )
567  {
568  // In the case that an initial condition is specified for an
569  // inductor, the DC op should be set up like a current source just
570  // for the operating point calculation.
571  dFdx[li_Pos][APosEquBraVarOffset] += 0.0;
572  dFdx[li_Neg][ANegEquBraVarOffset] += 0.0;
573  dFdx[li_Bra][ABraEquPosNodeOffset] += 0.0;
574  dFdx[li_Bra][ABraEquNegNodeOffset] += 0.0;
575  dFdx[li_Bra][ABraEquBraVarOffset] += 1.0;
576  }
577  else
578  {
579  dFdx[li_Pos][APosEquBraVarOffset] += 1.0;
580  dFdx[li_Neg][ANegEquBraVarOffset] -= 1.0;
581  dFdx[li_Bra][ABraEquPosNodeOffset] -= 1.0;
582  dFdx[li_Bra][ABraEquNegNodeOffset] += 1.0;
583  dFdx[li_Bra][ABraEquBraVarOffset] += 0.0;
584  }
585 
586  return true;
587 }
588 
589 //-----------------------------------------------------------------------------
590 // Function : Instance::setIC
591 // Purpose :
592 // Special Notes :
593 // Scope : public
594 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
595 // Creation Date : 1/11/02
596 //-----------------------------------------------------------------------------
598 {
599  int i_bra_sol;
600  int i_f_state;
601  double * nextStaVector = extData.nextStaVectorRawPtr;
602  double * currStaVector = extData.currStaVectorRawPtr;
603 
604  double * nextSolVector = extData.nextSolVectorRawPtr;
605  double * currSolVector = extData.currSolVectorRawPtr;
606 
607  if (ICGiven)
608  {
609  f0 = L*IC;
610  currStaVector[li_fstate] = f0;
611  nextStaVector[li_fstate] = f0;
612 
613  currSolVector[li_Bra] = IC;
614  nextSolVector[li_Bra] = IC;
615  }
616 
617  return true;
618 }
619 
620 //-----------------------------------------------------------------------------
621 // Function : Instance::varTypes
622 // Purpose :
623 // Special Notes :
624 // Scope : public
625 // Creator : Rob Hoekstra, SNL, Parallel Computational Sciences
626 // Creation Date : 2/17/04
627 //-----------------------------------------------------------------------------
628 void Instance::varTypes( std::vector<char> & varTypeVec )
629 {
630  varTypeVec.resize(1);
631  varTypeVec[0] = 'I';
632 }
633 
634 //-----------------------------------------------------------------------------
635 // Function : Model::processParams
636 // Purpose :
637 // Special Notes :
638 // Scope : public
639 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
640 // Creation Date : 6/03/02
641 //-----------------------------------------------------------------------------
643 {
644  return true;
645 }
646 
647 //----------------------------------------------------------------------------
648 // Function : Model::processInstanceParams
649 // Purpose :
650 // Special Notes :
651 // Scope : public
652 // Creator : Dave Shirely, PSSI
653 // Creation Date : 03/23/06
654 //----------------------------------------------------------------------------
656 {
657  std::vector<Instance*>::iterator iter;
658  std::vector<Instance*>::iterator first = instanceContainer.begin();
659  std::vector<Instance*>::iterator last = instanceContainer.end();
660 
661  for (iter=first; iter!=last; ++iter)
662  {
663  (*iter)->processParams();
664  }
665 
666  return true;
667 }
668 
669 //-----------------------------------------------------------------------------
670 // Function : Model::Model
671 // Purpose : block constructor
672 // Special Notes :
673 // Scope : public
674 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
675 // Creation Date : 3/16/00
676 //-----------------------------------------------------------------------------
678  const Configuration & configuration,
679  const ModelBlock & MB,
680  const FactoryBlock & factory_block)
681  : DeviceModel(MB, configuration.getModelParameters(), factory_block),
682  inductanceMultiplier(1.0),
683  IC(0.0),
684  tempCoeff1(0.0),
685  tempCoeff2(0.0),
686  tnom(getDeviceOptions().tnom),
687  tnomGiven(0)
688 {
689 
690  // Set params to constant default values:
691  setDefaultParams ();
692 
693  // Set params according to .model line and constant defaults from metadata:
694  setModParams (MB.params);
695 
696  // Set any non-constant parameter defaults:
697  if (!given("TNOM"))
699 
700  // Calculate any parameters specified as expressions:
702 
703  // calculate dependent (ie computed) params and check for errors:
704 
705  processParams ();
706 }
707 
708 //-----------------------------------------------------------------------------
709 // Function : Model::~Model
710 // Purpose : destructor
711 // Special Notes :
712 // Scope : public
713 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
714 // Creation Date : 3/16/00
715 //-----------------------------------------------------------------------------
717 {
718  std::vector<Instance*>::iterator iter;
719  std::vector<Instance*>::iterator first = instanceContainer.begin();
720  std::vector<Instance*>::iterator last = instanceContainer.end();
721 
722  for (iter=first; iter!=last; ++iter)
723  {
724  delete (*iter);
725  }
726 
727 }
728 
729 //-----------------------------------------------------------------------------
730 // Function : Model::printOutInstances
731 // Purpose : debugging tool.
732 // Special Notes :
733 // Scope : public
734 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
735 // Creation Date : 4/04/00
736 //-----------------------------------------------------------------------------
737 std::ostream &Model::printOutInstances(std::ostream &os) const
738 {
739  std::vector<Instance*>::const_iterator iter;
740  std::vector<Instance*>::const_iterator first = instanceContainer.begin();
741  std::vector<Instance*>::const_iterator last = instanceContainer.end();
742 
743  int i, isize;
744  isize = instanceContainer.size();
745 
746  os << std::endl;
747  os << "Number of Inductor instances: " << isize << std::endl;
748  os << " name=\t\tmodelName\tParameters" << std::endl;
749  for (i=0, iter=first; iter!=last; ++iter, ++i)
750  {
751  os << " " << i << ": " << (*iter)->getName() << "\t";
752  os << getName();
753  os << "\t\tL = " << (*iter)->L;
754  os << "\tIC = " << (*iter)->IC;
755  os << std::endl;
756  }
757 
758  os << std::endl;
759 
760  return os;
761 }
762 
763 //-----------------------------------------------------------------------------
764 // Function : Model::forEachInstance
765 // Purpose :
766 // Special Notes :
767 // Scope : public
768 // Creator : David Baur
769 // Creation Date : 2/4/2014
770 //-----------------------------------------------------------------------------
771 /// Apply a device instance "op" to all instances associated with this
772 /// model
773 ///
774 /// @param[in] op Operator to apply to all instances.
775 ///
776 ///
777 void Model::forEachInstance(DeviceInstanceOp &op) const /* override */
778 {
779  for (std::vector<Instance *>::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
780  op(*it);
781 }
782 
783 
784 // Inductor Master functions:
785 
786 //-----------------------------------------------------------------------------
787 // Function : Master::updateState
788 // Purpose :
789 // Special Notes :
790 // Scope : public
791 // Creator : Eric Keiter, SNL
792 // Creation Date : 11/26/08
793 //-----------------------------------------------------------------------------
794 bool Master::updateState (double * solVec, double * staVec, double * stoVec)
795 {
796  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
797  {
798  Instance & inst = *(*it);
799 
800  double current = solVec[inst.li_Bra];
801 
802  if( (getSolverState().dcopFlag) && inst.ICGiven )
803  current = inst.IC;
804 
805  inst.f0 = inst.L*current;
806  staVec[inst.li_fstate] = inst.f0;
807  }
808 
809  return true;
810 }
811 
812 //-----------------------------------------------------------------------------
813 // Function : Master::updateSecondaryState
814 // Purpose :
815 // Special Notes :
816 // Scope : public
817 // Creator : Eric Keiter, SNL
818 // Creation Date : 11/26/08
819 //-----------------------------------------------------------------------------
820 bool Master::updateSecondaryState ( double * staDerivVec, double * stoVec )
821 {
822  return true;
823 }
824 
825 //-----------------------------------------------------------------------------
826 // Function : Master::loadDAEVectors
827 // Purpose :
828 // Special Notes :
829 // Scope : public
830 // Creator : Eric Keiter, SNL
831 // Creation Date : 11/26/08
832 //-----------------------------------------------------------------------------
833 bool Master::loadDAEVectors (double * solVec, double * fVec, double *qVec, double * bVec, double * storeLeadF, double * storeLeadQ)
834 {
835  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
836  {
837  Instance & inst = *(*it);
838 
839  double current = 0.0;
840  double coef = 0.0;
841 
842  double v_pos = solVec[inst.li_Pos];
843  double v_neg = solVec[inst.li_Neg];
844  double vind = v_pos-v_neg;
845 
846  // In the case that an initial condition is specified, the inductor is set up
847  // like a current source for the DC operating point. We don't deal with the
848  // node voltages in that case, so set coef to 0.
849  if (getSolverState().dcopFlag && inst.ICGiven)
850  {
851  current = inst.IC;
852  solVec[inst.li_Bra] = current;
853  coef = 0.0;
854  }
855  else
856  {
857  current = solVec[inst.li_Bra];
858  coef = -vind;
859  }
860 
861  // load the current into the two KCL rhs vector elements
862  fVec[inst.li_Pos] += current;
863  fVec[inst.li_Neg] += -current;
864  fVec[inst.li_Bra] += coef;
865  qVec[inst.li_Bra] += inst.f0;
866  }
867 return true;
868 }
869 
870 //-----------------------------------------------------------------------------
871 // Function : Master::loadDAEMatrices
872 // Purpose :
873 // Special Notes :
874 // Scope : public
875 // Creator : Eric Keiter, SNL
876 // Creation Date : 11/26/08
877 //-----------------------------------------------------------------------------
878 bool Master::loadDAEMatrices (N_LAS_Matrix & dFdx, N_LAS_Matrix & dQdx)
879 {
880  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
881  {
882  Instance & inst = *(*it);
883 
884  if ( getSolverState().dcopFlag && inst.ICGiven )
885  {
886  // In the case that an initial condition is specified for an
887  // inductor, the DC op should be set up like a current source just
888  // for the operating point calculation.
889 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
890  *inst.fBraEquBraVarPtr += 1.0;
891 #else
892  dFdx[inst.li_Bra][inst.ABraEquBraVarOffset] += 1.0;
893 #endif
894  }
895  else
896  {
897 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
898  *inst.fPosEquBraVarPtr += 1.0;
899  *inst.fNegEquBraVarPtr -= 1.0;
900  *inst.fBraEquPosNodePtr -= 1.0;
901  *inst.fBraEquNegNodePtr += 1.0;
902 #else
903  dFdx[inst.li_Pos][inst.APosEquBraVarOffset] += 1.0;
904  dFdx[inst.li_Neg][inst.ANegEquBraVarOffset] -= 1.0;
905  dFdx[inst.li_Bra][inst.ABraEquPosNodeOffset] -= 1.0;
906  dFdx[inst.li_Bra][inst.ABraEquNegNodeOffset] += 1.0;
907 #endif
908  }
909 
910 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
911  *inst.qBraEquBraVarPtr += inst.L;
912 #else
913  dQdx[inst.li_Bra][inst.ABraEquBraVarOffset] += inst.L;
914 #endif
915  }
916  return true;
917 }
918 
919 Device *Traits::factory(const Configuration &configuration, const FactoryBlock &factory_block)
920 {
921 
922  return new Master(configuration, factory_block, factory_block.solverState_, factory_block.deviceOptions_);
923 }
924 
926 {
928  .registerDevice("l", 1)
929  .registerModelType("l", 1)
930  .registerModelType("ind", 1);
931 }
932 
933 //-----------------------------------------------------------------------------
934 // Function : indSensitivity::operator
935 // Purpose : produces df/dp and dq/dp, where p=L.
936 // Special Notes :
937 // Scope : public
938 // Creator : Eric Keiter, SNL
939 // Creation Date : 7/31/2014
940 //-----------------------------------------------------------------------------
942  const ParameterBase &entity,
943  const std::string & name,
944  std::vector<double> & dfdp,
945  std::vector<double> & dqdp,
946  std::vector<double> & dbdp,
947  std::vector<int> & Findices,
948  std::vector<int> & Qindices,
949  std::vector<int> & Bindices
950  ) const
951 {
952  const ParameterBase * e1 = &entity;
953  const Instance * in = dynamic_cast<const Instance *> (e1);
954 
955  double * solVec = in->extData.nextSolVectorRawPtr;
956  double current = solVec[in->li_Bra];
957  if( (in->getSolverState().dcopFlag) && in->ICGiven )
958  {
959  current = in->IC;
960  }
961 
962  double dqdpLoc = current;
963 
964  dqdp.resize(1);
965  dqdp[0] = dqdpLoc;
966 
967  Qindices.resize(1);
968  Qindices[0] = in->li_Bra;
969 }
970 
971 } // namespace Inductor
972 } // namespace Device
973 } // namespace Xyce