Xyce  6.1
N_LOA_HBLoader.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_LOA_HBLoader.C,v $
27 // Purpose :
28 // Special Notes :
29 // Creator :
30 // Creation Date :
31 //
32 // Revision Information:
33 // ---------------------
34 // Revision Number: $Revision: 1.48.2.1 $
35 // Revision Date : $Date: 2015/04/02 18:20:16 $
36 // Current Owner : $Author: tvrusso $
37 //-----------------------------------------------------------------------------
38 
39 #include <Xyce_config.h>
40 
41 
42 // ---------- Standard Includes ----------
43 #include <iostream>
44 
45 // ---------- Xyce Includes ----------
46 
47 #include <N_LOA_HBLoader.h>
48 #include <N_MPDE_Discretization.h>
49 #include <N_ERH_ErrorMgr.h>
50 
51 #include <N_LAS_Builder.h>
52 #include <N_LAS_HBBuilder.h>
53 #include <N_LAS_Matrix.h>
54 #include <N_LAS_BlockVector.h>
55 #include <N_LAS_BlockMatrix.h>
56 #include <N_LAS_BlockSystemHelpers.h>
57 
58 #include <N_UTL_fwd.h>
59 #include <N_UTL_FeatureTest.h>
60 #include <N_UTL_Math.h>
61 
62 #include <N_PDS_ParMap.h>
63 #include <N_DEV_DeviceMgr.h>
64 
65 #include <Epetra_MultiVector.h>
66 #include <Epetra_BlockMap.h>
67 
68 using Teuchos::rcp;
69 using Teuchos::RCP;
70 using Teuchos::rcp_dynamic_cast;
71 
72 namespace Xyce {
73 namespace Loader {
74 
75 
76  // Default constructor
78  const Teuchos::RCP<const N_MPDE_Discretization> discPtr,
79  Device::DeviceMgr & device_manager,
80  Linear::Builder & builder)
81  : CktLoader(device_manager),
82  deviceManager_(device_manager),
83  builder_(builder),
84  fastTimeDiscPtr_(discPtr),
85  periodicTimesOffset_(0),
86  period_(0),
87  matrixFreeFlag_(false)
88 {
89  // Now initialize all the time domain working vectors.
90  appVecPtr_ = rcp(builder_.createVector());
91  appNextStaVecPtr_ = rcp(builder_.createStateVector());
92  appCurrStaVecPtr_ = rcp(builder_.createStateVector());
93  appLastStaVecPtr_ = rcp(builder_.createStateVector());
94 
95  appdQdxPtr_ = rcp(builder_.createMatrix());
96  appdFdxPtr_ = rcp(builder_.createMatrix());
97 
98  appNextStoVecPtr_ = rcp(builder_.createStoreVector());
99  appCurrStoVecPtr_ = rcp(builder_.createStoreVector());
100  appLastStoVecPtr_ = rcp(builder_.createStoreVector());
101  appStoLeadCurrQVecPtr_ = rcp(builder_.createStoreVector());
102 
103  appNextLeadFVecPtr_ = rcp(builder.createLeadCurrentVector());
104  appCurrLeadFVecPtr_ = rcp(builder.createLeadCurrentVector());
105  appLastLeadFVecPtr_ = rcp(builder.createLeadCurrentVector());
106  appLeadQVecPtr_ = rcp(builder.createLeadCurrentVector());
107  appNextJunctionVVecPtr_ = rcp(builder.createLeadCurrentVector());
108  appCurrJunctionVVecPtr_ = rcp(builder.createLeadCurrentVector());
109  appLastJunctionVVecPtr_ = rcp(builder.createLeadCurrentVector());
110 }
111 
112 //-----------------------------------------------------------------------------
113 // Function : HBLoader::registerHBBuilder
114 // Purpose : Registration method for the HB builder
115 // Special Notes :
116 // Scope : public
117 // Creator : Heidi Thornquist, Sandia Labs
118 // Creation Date : 11/08/13
119 //-----------------------------------------------------------------------------
120 void HBLoader::registerHBBuilder( Teuchos::RCP<Linear::HBBuilder> hbBuilderPtr )
121 {
122  hbBuilderPtr_ = hbBuilderPtr;
123 
124  // Now initialize all the frequency domain working vectors.
125  bXtPtr_ = hbBuilderPtr_->createTimeDomainBlockVector();
126  bVtPtr_ = hbBuilderPtr_->createTimeDomainBlockVector();
127  bmdQdxPtr_ = rcp_dynamic_cast<Linear::BlockMatrix>(rcp(hbBuilderPtr_->createMatrix()));
128  bmdFdxPtr_ = rcp_dynamic_cast<Linear::BlockMatrix>(rcp(hbBuilderPtr_->createMatrix()));
129 
130  // Vectors related to lead currents
131  bStoreVecFreqPtr_ = hbBuilderPtr_->createExpandedRealFormTransposeStoreBlockVector();
132  bStoreLeadCurrQVecFreqPtr_ = hbBuilderPtr_->createExpandedRealFormTransposeStoreBlockVector();
133 
134  bLeadCurrentVecFreqPtr_ = hbBuilderPtr_->createExpandedRealFormTransposeLeadCurrentBlockVector();
135  bLeadCurrentQVecFreqPtr_ = hbBuilderPtr_->createExpandedRealFormTransposeLeadCurrentBlockVector();
136 }
137 
138 //-----------------------------------------------------------------------------
139 // Function : HBLoader::setFreq
140 // Purpose : Assign times for fast time scale
141 // Special Notes :
142 // Scope : public
143 // Creator :
144 // Creation Date :
145 //-----------------------------------------------------------------------------
146 void HBLoader::setHBFreqs( const std::vector<double> & freqs)
147 {
148  freqs_ = freqs;
149 }
150 
151 
152 
153 //-----------------------------------------------------------------------------
154 // Function : HBLoader::setFastTimes
155 // Purpose : Assign times for fast time scale
156 // Special Notes :
157 // Scope : public
158 // Creator :
159 // Creation Date :
160 //-----------------------------------------------------------------------------
161 void HBLoader::setFastTimes( const std::vector<double> & times )
162 {
163  times_ = times;
164 }
165 
166 //-----------------------------------------------------------------------------
167 // Function : HBLoader::loadDAEMatrices
168 // Purpose :
169 // Special Notes :
170 // Scope : public
171 // Creator :
172 // Creation Date :
173 //-----------------------------------------------------------------------------
174 bool HBLoader::loadDAEMatrices( Linear::Vector * X,
175  Linear::Vector * S,
176  Linear::Vector * dSdt,
177  Linear::Vector * Store,
178  Linear::Matrix * dQdx,
179  Linear::Matrix * dFdx)
180 {
181  if ( matrixFreeFlag_ )
182  {
183  if (DEBUG_HB)
184  {
185  Xyce::dout() << std::endl
186  << Xyce::section_divider << std::endl
187  << " HBLoader::loadDAEMatrices: matrixFree case" << std::endl;
188  }
189 
190  dQdx->put(0.0);
191  dFdx->put(0.0);
192  // Do nothing in the Matrix Free Case.
193  return(true);
194  }
195  else
196  {
197  if (DEBUG_HB)
198  {
199  Xyce::dout() << std::endl
200  << Xyce::section_divider << std::endl
201  << " HBLoader::loadDAEMatrices: Time dependent with matrix case" << std::endl;
202  }
203 
204  return loadTimeDepDAEMatrices(X,S,dSdt, Store, dQdx,dFdx);
205  }
206 
207 }
208 
209 //-----------------------------------------------------------------------------
210 // Function : HBLoader::loadTimeDepDAEMatrices
211 // Purpose :
212 // Special Notes :
213 // Scope : public
214 // Creator : Todd Coffey, 1414
215 // Creation Date : 09/03/08
216 //-----------------------------------------------------------------------------
217 bool HBLoader::loadTimeDepDAEMatrices( Linear::Vector * X,
218  Linear::Vector * S,
219  Linear::Vector * dSdt,
220  Linear::Vector * Store,
221  Linear::Matrix * dQdx,
222  Linear::Matrix * dFdx)
223 {
224  if (DEBUG_HB)
225  {
226  Xyce::dout() << std::endl
227  << Xyce::section_divider << std::endl
228  << " HBLoader::loadTimeDepDAEMatrices" << std::endl;
229  }
230 
231  //Zero out matrices
232  dQdx->put(0.0);
233  dFdx->put(0.0);
234 
235  Linear::Vector appdSdt( *appNextStaVecPtr_ );
236 
237  Linear::BlockMatrix & bdQdx = *dynamic_cast<Linear::BlockMatrix*>(dQdx);
238  Linear::BlockMatrix & bdFdx = *dynamic_cast<Linear::BlockMatrix*>(dFdx);
239  Linear::BlockVector & bX = *dynamic_cast<Linear::BlockVector*>(X);
240 #ifdef Xyce_FLEXIBLE_DAE_LOADS
241  Linear::BlockVector & bS = *dynamic_cast<Linear::BlockVector*>(S);
242  Linear::BlockVector & bdSdt = *dynamic_cast<Linear::BlockVector*>(dSdt);
243  Linear::BlockVector & bStore = *dynamic_cast<Linear::BlockVector*>(Store);
244 #endif // Xyce_FLEXIBLE_DAE_LOADS
245 
246  int BlockCount = bX.blockCount();
247 
248  for( int i = 0; i < BlockCount; ++i )
249  {
250  if (DEBUG_HB)
251  {
252  Xyce::dout() << "Processing diagonal matrix block " << i << " of " << BlockCount-1 << std::endl;
253  }
254 
255 #ifdef Xyce_FLEXIBLE_DAE_LOADS
256  //Set Time for fast time scale somewhere
257 // state_.fastTime = times_[i];
259 
260  //Update the sources
261  appLoaderPtr_->updateSources();
262 
263  *appVecPtr_ = bX.block(i);
264  *appNextStaVecPtr_ = bS.block(i);
265  appdSdt = bdSdt.block(i);
266  *appNextStoVecPtr_ = bStore.block(i);
267 
268  appLoaderPtr_->loadDAEMatrices( &*appVecPtr_, &*appNextStaVecPtr_, &appdSdt, &*appNextStoVecPtr_, &*appdQdxPtr_, &*appdFdxPtr_);
269 
270  bdQdx.block(i,i).add( *appdQdxPtr_ );
271  bdFdx.block(i,i).add( *appdFdxPtr_ );
272 #else
273  //For now, the matrices are loaded during the loadDAEVectors method
274  //Just copied here
275  bdQdx.block(i,i).add( bmdQdxPtr_->block(i,i) );
276  bdFdx.block(i,i).add( bmdFdxPtr_->block(i,i) );
277 
278 #endif // Xyce_FLEXIBLE_DAE_LOADS
279  }
280 
281  // Fast Time scale discretization terms:
282  // These are d(dQ/dt1)/dx terms, but go into bdFdx.
283  // For this procedure, need to re-use the app matrix, appdQdx.
284  Linear::Matrix & dQdxdt = *appdQdxPtr_;
285 
286  int Start = fastTimeDiscPtr_->Start();
287  int Width = fastTimeDiscPtr_->Width();
288 
289  const std::vector<double> & Coeffs = fastTimeDiscPtr_->Coeffs();
290 
291  for( int i = 0; i < BlockCount; ++i )
292  {
293  if (DEBUG_HB)
294  {
295  Xyce::dout() << "Processing off diagonal matrix blocks on row " << i << " of " << BlockCount-1 << std::endl;
296  }
297 
298  int Loc;
299  int indexT1 = i + Start + periodicTimesOffset_;
300  int indexT2 = indexT1 + Width - 1;
301  double invh2 = 1.0 / (periodicTimes_[indexT2] - periodicTimes_[indexT1]);
302 
303  for( int j = 0; j < Width; ++j )
304  {
305  Loc = i + (j+Start);
306 
307  if( Loc < 0 )
308  {
309  Loc += BlockCount;
310  }
311  else if( Loc > (BlockCount-1) )
312  {
313  Loc -= BlockCount;
314  }
315 
316  dQdxdt.put(0.0);
317  dQdxdt.add( bdQdx.block(Loc,Loc) );
318  dQdxdt.scale( Coeffs[j]*invh2 );
319  bdFdx.block(i,Loc).add( dQdxdt );
320  }
321  }
322 
323  dQdx->fillComplete();
324  dFdx->fillComplete();
325 
326  if (DEBUG_HB)
327  {
328  Xyce::dout() << "HB bX:" << std::endl;
329  bX.printPetraObject(dout());
330  Xyce::dout() << "HB bdQdx:" << std::endl;
331  bdQdx.printPetraObject(dout());
332  Xyce::dout() << "HB bdFdx:" << std::endl;
333  bdFdx.printPetraObject(dout());
334 #ifdef Xyce_FLEXIBLE_DAE_LOADS
335  Xyce::dout() << "HB bS:" << std::endl;
336  bS.printPetraObject(dout());
337  Xyce::dout() << "HB dSdt:" << std::endl;
338  bdSdt.printPetraObject(dout());
339  Xyce::dout() << "HB bStore:" << std::endl;
340  bStore.printPetraObject(dout());
341 #endif // Xyce_FLEXIBLE_DAE_LOADS
342 
343  Xyce::dout() << Xyce::section_divider << std::endl;
344  }
345 
346  return true;
347 }
348 
349 //-----------------------------------------------------------------------------
350 // Function : HBLoader::applyDAEMatrices
351 // Purpose :
352 // Special Notes :
353 // Scope : public
354 // Creator :
355 // Creation Date :
356 //-----------------------------------------------------------------------------
357 bool HBLoader::applyDAEMatrices( Linear::Vector * Xf,
358  Linear::Vector * S,
359  Linear::Vector * dSdt,
360  Linear::Vector * Store,
361  const Linear::Vector & Vf,
362  Linear::Vector * dQdxV,
363  Linear::Vector * dFdxV )
364 {
365  if ( !matrixFreeFlag_ )
366  {
367  std::string msg="HBLoader::applyDAEMatrices. This function should only be called in the matrix free case.";
368  N_ERH_ErrorMgr::report ( N_ERH_ErrorMgr::DEV_FATAL, msg);
369  }
370  if (DEBUG_HB)
371  {
372  Xyce::dout() << std::endl
373  << Xyce::section_divider << std::endl
374  << " HBLoader::applyDAEMatrices" << std::endl;
375  }
376 
377  //Zero out matvec vectors
378  dQdxV->putScalar(0.0); // This is dQdx * V on output (only used in HB-env simulation)
379  dFdxV->putScalar(0.0); // This is dFdx * V on output (This is main output for HB)
380 
381  bXtPtr_->putScalar(0.0);
382  bVtPtr_->putScalar(0.0);
383 
384  Linear::BlockVector & bXf = *dynamic_cast<Linear::BlockVector*>(Xf);
385  // We have to do something special with Vf because AztecOO (or Belos)
386  // probably used the Epetra_LinearProblem's Epetra_Maps to create the input
387  // vector here. In this case, Vf is just an Linear::Vector and not an
388  // Linear::BlockVector.
389  const Linear::BlockVector bVf(Vf, bXf.blockSize());
390 
391  permutedIFT(bXf, &*bXtPtr_);
392  permutedIFT(bVf, &*bVtPtr_);
393 
394  Linear::BlockVector & bX = *bXtPtr_;
395 
396  Linear::BlockVector * bdQdxV = dynamic_cast<Linear::BlockVector*>(dQdxV);
397  Linear::BlockVector * bdFdxV = dynamic_cast<Linear::BlockVector*>(dFdxV);
398 
399  Teuchos::RCP<Linear::BlockVector> bdQdxVt = hbBuilderPtr_->createTimeDomainBlockVector();
400  Teuchos::RCP<Linear::BlockVector> bdFdxVt = hbBuilderPtr_->createTimeDomainBlockVector();
401 
402  int BlockCount = bX.blockCount();
403  for( int i = 0; i < BlockCount; ++i )
404  {
405  if (DEBUG_HB)
406  {
407  Xyce::dout() << "Processing diagonal matrix block " << i << " of " << BlockCount-1 << std::endl;
408  }
409 
410  // Get the already stored time-domain Jacobian matrices
413 
414  if (DEBUG_HB)
415  {
416  Xyce::dout() << "bVtPtr block i = " << i << " : " << std::endl;
417  bVtPtr_->block(i).printPetraObject(dout());
418 
419  Xyce::dout() << "appdQdxPtr_ = " << i << " : " << std::endl;
420  appdQdxPtr_->printPetraObject(dout());
421 
422  Xyce::dout() << "appdFdxPtr_ = " << i << " : " << std::endl;
423  appdFdxPtr_->printPetraObject(dout());
424  }
425 
426 
427  appdQdxPtr_->matvec(false, bVtPtr_->block(i), bdQdxVt->block(i));
428  appdFdxPtr_->matvec(false, bVtPtr_->block(i), bdFdxVt->block(i));
429 
430  if (DEBUG_HB)
431  {
432  Xyce::dout() << "bdQdxVt block i = " << i << " : " << std::endl;
433  bdQdxVt->block(i).printPetraObject(dout());
434 
435  Xyce::dout() << "bdFdxVt block i = " << i << " : " << std::endl;
436  bdFdxVt->block(i).printPetraObject(dout());
437  }
438  }
439 
440  permutedFFT(*bdQdxVt, bdQdxV);
441  permutedFFT(*bdFdxVt, bdFdxV);
442 
443  int blockCount = bXf.blockCount();
444  int blockSize = bXf.block(0).globalLength();
445 // double omega = 2.0 * M_PI/ period_;
446 
447  int size_ = freqs_.size();
448  int posFreq = (size_-1)/2;
449  double omega = 2.0 * M_PI * freqs_[posFreq];
450  for( int i = 0; i < blockCount; ++i )
451  {
452  Linear::Vector QVec(bXf.block(i));
453  Linear::Vector freqVec = bdQdxV->block(i);
454 
455  // Only one processor owns each block of the frequency-domain vector
456  if (freqVec.localLength() > 0)
457  {
458 
459  omega = 2.0 * M_PI * freqs_[posFreq];
460 
461  QVec[0] = -freqVec[1]*omega;
462  QVec[1] = freqVec[0]*omega;
463 
464  for (int j=1; j < (blockSize/2+1)/2; ++j)
465  {
466  omega = 2.0 * M_PI * freqs_[posFreq+j];
467  QVec[2*j] = -freqVec[2*j+1]*omega;
468  QVec[2*(blockSize/2-j)] = -freqVec[2*j+1]*omega;
469 
470  QVec[2*j+1] = freqVec[2*j]*omega;
471  QVec[2*(blockSize/2-j)+1] = -freqVec[2*j]*omega;
472  }
473  }
474 
475  bdFdxV->block(i).update(1.0, QVec , 1.0);
476  }
477 
478  if (DEBUG_HB)
479  {
480  Xyce::dout() << "HB bX:" << std::endl;
481  bX.printPetraObject(dout());
482  Xyce::dout() << "HB bdQdxV:" << std::endl;
483  bdQdxV->printPetraObject(dout());
484  Xyce::dout() << "HB bdFdxV:" << std::endl;
485  bdFdxV->printPetraObject(dout());
486 
487  Xyce::dout() << Xyce::section_divider << std::endl;
488  }
489 
490  return true;
491 }
492 
493 //-----------------------------------------------------------------------------
494 // Function : HBLoader::updateState
495 // Purpose :
496 // Special Notes : ERK. This function needs to be a no-op. The reason
497 // is that the state information needs to be the same
498 // at the time of updateState, loadDAEVectors and
499 // loadDAEMatrices. Thus, they have to all happen inside
500 // of the same "fast time" loop. So, this functionality
501 // has been moved up into loadDAEVectors.
502 //
503 // Note: for similar reasons, loadDAEMatrices is called from
504 // within that function as well.
505 //
506 // Scope : public
507 // Creator : Todd Coffey, 1414
508 // Creation Date : 01/17/07
509 //-----------------------------------------------------------------------------
511  (Linear::Vector * nextSolVectorPtr,
512  Linear::Vector * currSolVectorPtr,
513  Linear::Vector * lastSolVectorPtr,
514  Linear::Vector * nextStaVectorPtr,
515  Linear::Vector * currStaVectorPtr,
516  Linear::Vector * lastStaVectorPtr,
517  Linear::Vector * nextStoVectorPtr,
518  Linear::Vector * currStoVectorPtr,
519  Linear::Vector * lastStoVectorPtr
520  )
521 {
522  bool bsuccess = true;
523 
524  // For HB case, this needs to be a no-op.
525 
526  return bsuccess;
527 }
528 
529 //-----------------------------------------------------------------------------
530 // Function : HBLoader::loadDAEVectors
531 // Purpose :
532 // Special Notes :
533 // Scope : public
534 // Creator :
535 // Creation Date :
536 //-----------------------------------------------------------------------------
537 bool HBLoader::loadDAEVectors( Linear::Vector * Xf,
538  Linear::Vector * currX,
539  Linear::Vector * lastX,
540  Linear::Vector * S,
541  Linear::Vector * currS,
542  Linear::Vector * lastS,
543  Linear::Vector * dSdt,
544  Linear::Vector * Store,
545  Linear::Vector * currStore,
546  Linear::Vector * lastStore,
547  Linear::Vector * storeLeadCurrQ,
548  Linear::Vector * nextLeadFVectorPtr,
549  Linear::Vector * currLeadFVectorPtr,
550  Linear::Vector * lastLeadFVectorPtr,
551  Linear::Vector * nextLeadQVectorPtr,
552  Linear::Vector * nextJunctionVVectorPtr,
553  Linear::Vector * currJunctionVVectorPtr,
554  Linear::Vector * lastJunctionVVectorPtr,
555  Linear::Vector * Q,
556  Linear::Vector * F,
557  Linear::Vector * B,
558  Linear::Vector * dFdxdVp,
559  Linear::Vector * dQdxdVp )
560 {
561  if (DEBUG_HB)
562  {
563  Xyce::dout() << std::endl
564  << Xyce::section_divider << std::endl
565  << " HBLoader::loadDAEVectors" << std::endl;
566  }
567 
568  //Zero out vectors
569  appVecPtr_->putScalar(0.0);
570  appNextStaVecPtr_->putScalar(0.0);
571  appCurrStaVecPtr_->putScalar(0.0);
572  appLastStaVecPtr_->putScalar(0.0);
573  Linear::Vector appdSdt( *appNextStaVecPtr_ );
574 
575  appNextStoVecPtr_->putScalar(0.0);
576  appCurrStoVecPtr_->putScalar(0.0);
577  appLastStoVecPtr_->putScalar(0.0);
578  appStoLeadCurrQVecPtr_->putScalar(0.0);
579 
580  Linear::Vector appQ( *appVecPtr_ );
581  Linear::Vector appF( *appVecPtr_ );
582  Linear::Vector appB( *appVecPtr_ );
583 
584  Linear::Vector appdFdxdVp( *appVecPtr_ );
585  Linear::Vector appdQdxdVp( *appVecPtr_ );
586 
587  // This is a temporary load storage vector.
588  Linear::Vector dQdt2( *appVecPtr_ );
589 
590  bXtPtr_->putScalar(0.0);
591 
592  Linear::BlockVector & bXf = *dynamic_cast<Linear::BlockVector*>(Xf);
593 
594  permutedIFT(bXf, &*bXtPtr_);
595 
596  // 12/8/06 tscoffe: Note: "b" at beginning of variable name means Linear::BlockVector
597  Linear::BlockVector & bX = *bXtPtr_;
598  Linear::BlockVector & bS = *dynamic_cast<Linear::BlockVector*>(S);
599  Linear::BlockVector & bcurrS = *dynamic_cast<Linear::BlockVector*>(currS);
600  Linear::BlockVector & blastS = *dynamic_cast<Linear::BlockVector*>(lastS);
601  Linear::BlockVector & bdSdt = *dynamic_cast<Linear::BlockVector*>(dSdt);
602  Linear::BlockVector & bStore = *dynamic_cast<Linear::BlockVector*>(Store);
603  Linear::BlockVector & bcurrStore = *dynamic_cast<Linear::BlockVector*>(currStore);
604  Linear::BlockVector & blastStore = *dynamic_cast<Linear::BlockVector*>(lastStore);
605 
606  Linear::BlockVector & bNextLeadF = *dynamic_cast<Linear::BlockVector*>(nextLeadFVectorPtr);
607  Linear::BlockVector & bCurrLeadF = *dynamic_cast<Linear::BlockVector*>(currLeadFVectorPtr);
608  Linear::BlockVector & bLastLeadF = *dynamic_cast<Linear::BlockVector*>(lastLeadFVectorPtr);
609  Linear::BlockVector & bLeadQ = *dynamic_cast<Linear::BlockVector*>(nextLeadQVectorPtr);
610  Linear::BlockVector & bNextJunctionV = *dynamic_cast<Linear::BlockVector*>(nextJunctionVVectorPtr);
611  Linear::BlockVector & bCurrJunctionV = *dynamic_cast<Linear::BlockVector*>(currJunctionVVectorPtr);
612  Linear::BlockVector & bLastJunctionV = *dynamic_cast<Linear::BlockVector*>(lastJunctionVVectorPtr);
613 
614  Linear::BlockVector * bQ = dynamic_cast<Linear::BlockVector*>(Q);
615  Linear::BlockVector * bF = dynamic_cast<Linear::BlockVector*>(F);
616  Linear::BlockVector * bB = dynamic_cast<Linear::BlockVector*>(B);
617 
618  Linear::BlockVector * bdFdxdVp = dynamic_cast<Linear::BlockVector*>(dFdxdVp);
619  Linear::BlockVector * bdQdxdVp = dynamic_cast<Linear::BlockVector*>(dQdxdVp);
620 
621 
622 // Linear::Vector * storeLeadCurrQ,
623 
624  Linear::BlockVector & bstoreLeadCurrQ = *dynamic_cast<Linear::BlockVector*>(storeLeadCurrQ);
625 
626  Linear::BlockVector bQt(*bXtPtr_);
627  Linear::BlockVector bFt(*bXtPtr_);
628  Linear::BlockVector bBt(*bXtPtr_);
629 
630  Linear::BlockVector bdQdxdVpt(*bXtPtr_);
631  Linear::BlockVector bdFdxdVpt(*bXtPtr_);
632 
633 
634 // Linear::BlockVector bstoreLeadCurrQt(*bXtPt );
635 
636 #ifndef Xyce_FLEXIBLE_DAE_LOADS
637  bmdQdxPtr_->put(0.0);
638  bmdFdxPtr_->put(0.0);
639 #endif
640 
641  int BlockCount = bX.blockCount();
642  int BlockSize = bX.blockSize();
643 
644 #if 0
645  Xyce::dout() << "bX.blockCount = "<< BlockCount <<std::endl;
646  Xyce::dout() << "bX.blockSize = "<< BlockSize <<std::endl;
647 #endif
648 
649  // We are storing the time domain Jacobians, initialize the memory here
650  if ((int)vecAppdQdxPtr_.size()!=BlockCount)
651  {
652  vecAppdQdxPtr_.resize( BlockCount );
653  vecAppdFdxPtr_.resize( BlockCount );
654  for ( int i=0; i < BlockCount; ++i )
655  {
656  vecAppdQdxPtr_[i] = rcp(builder_.createMatrix());
657  vecAppdFdxPtr_[i] = rcp(builder_.createMatrix());
658  }
659  }
660 
661  // Now perform implicit application of frequency domain Jacobian.
662  for( int i = 0; i < BlockCount; ++i )
663  {
664  if (DEBUG_HB)
665  {
666  Xyce::dout() << "Processing vectors for block " << i << " of " << BlockCount-1 << std::endl;
667  }
668 
669  //Set Time for fast time scale somewhere
670 // state_.fastTime = times_[i];
672 
673  if (DEBUG_HB)
674  {
675  Xyce::dout() << "Calling updateSources on the appLoader" << std::endl;
676  }
677 
678  //Update the sources
679  appLoaderPtr_->updateSources(); // this is here to handle "fast" sources.
680 
681  *appVecPtr_ = bX.block(i);
682  *appNextStaVecPtr_ = bS.block(i);
683  *appCurrStaVecPtr_ = bcurrS.block(i);
684  *appLastStaVecPtr_ = blastS.block(i);
685  appdSdt = bdSdt.block(i);
686  *appNextStoVecPtr_ = bStore.block(i);
687  *appCurrStoVecPtr_ = bcurrStore.block(i);
688  *appLastStoVecPtr_ = blastStore.block(i);
689  *appStoLeadCurrQVecPtr_ = bstoreLeadCurrQ.block(i);
690  *appNextLeadFVecPtr_ = bNextLeadF.block(i);
691  *appCurrLeadFVecPtr_ = bCurrLeadF.block(i);
692  *appLastLeadFVecPtr_ = bLastLeadF.block(i);
693  *appLeadQVecPtr_ = bLeadQ.block(i);
694  *appNextJunctionVVecPtr_ = bNextJunctionV.block(i);
695  *appCurrJunctionVVecPtr_ = bCurrJunctionV.block(i);
696  *appLastJunctionVVecPtr_ = bLastJunctionV.block(i);
697 
698 //blastStore.block(i); // Need to update to correct block!
699 
700  if (DEBUG_HB)
701  {
702  Xyce::dout() << "Updating State for block " << i << " of " << BlockCount-1 << std::endl;
703  }
704 
705  // Note: This updateState call is here (instead of in the
706  // HBLoader::updateState function) because it has to be called
707  // for the same fast time point.
708  appLoaderPtr_->updateState
709  ( &*appVecPtr_,
710  &*appVecPtr_, // note, this is a placeholder! ERK
711  &*appVecPtr_, // note, this is a placeholder! ERK
714  );
715 
716  bS.block(i) = *appNextStaVecPtr_;
717  bcurrS.block(i) = *appCurrStaVecPtr_;
718  blastS.block(i) = *appLastStaVecPtr_;
719  bStore.block(i) = *appNextStoVecPtr_;
720  bcurrStore.block(i) = *appCurrStoVecPtr_;
721  blastStore.block(i) = *appLastStoVecPtr_;
722 
723  if (DEBUG_HB)
724  {
725  Xyce::dout() << "Calling loadDAEVectors on the appLoader" << std::endl;
726  }
727 
728  // This has to be done because the app loader does NOT zero these vectors out.
729  appQ.putScalar(0.0);
730  appF.putScalar(0.0);
731  appB.putScalar(0.0);
732  appdFdxdVp.putScalar(0.0);
733  appdQdxdVp.putScalar(0.0);
734 
735  appLoaderPtr_->loadDAEVectors
736  ( &*appVecPtr_,
737  &*appVecPtr_, // note, this is a placeholder! ERK
738  &*appVecPtr_, // note, this is a placeholder! ERK
743  &appQ, &appF, &appB,
744  &appdFdxdVp, &appdQdxdVp );
745 
746  bQt.block(i) = appQ;
747  bFt.block(i) = appF;
748  bBt.block(i) = appB;
749 
750  bdQdxdVpt.block(i) = appdQdxdVp;
751  bdFdxdVpt.block(i) = appdFdxdVp;
752 
753  bstoreLeadCurrQ.block(i) = *appStoLeadCurrQVecPtr_;
754 
755  bStore.block(i) = *appNextStoVecPtr_; //lead current get loaded during loadDAEVectors
756 
757  bNextLeadF.block(i) = *appNextLeadFVecPtr_; // lead currents loaded into lead current vector.
758  bLeadQ.block(i) = *appLeadQVecPtr_;
759 
760 // Xyce::dout() << "HB Loader Store Vector TD" << std::endl;
761 // bStore.printPetraObject(std::cout);
762 
763 // Xyce::dout() << "HB Loader Store Q Vector TD" << std::endl;
764 // bstoreLeadCurrQ.printPetraObject(std::cout);
765 
766  // Store the time domain Jacobian for future use.
767  vecAppdQdxPtr_[i]->put(0.0);
768  vecAppdFdxPtr_[i]->put(0.0);
769 
770  // Load dQdx and dFdx into the storage location for this time point.
771  appLoaderPtr_->loadDAEMatrices( &*appVecPtr_, &*appNextStaVecPtr_, &appdSdt, &*appNextStoVecPtr_, &*vecAppdQdxPtr_[i], &*vecAppdFdxPtr_[i]);
772  }
773 
774  permutedFFT(bBt, bB);
775 
776  permutedFFT(bQt, bQ);
777 
778  permutedFFT(bFt, bF);
779 
780  permutedFFT(bdQdxdVpt, bdQdxdVp);
781 
782  permutedFFT(bdFdxdVpt, bdFdxdVp);
783 
784  bStoreVecFreqPtr_->putScalar(0.0);
785  bStoreLeadCurrQVecFreqPtr_->putScalar(0.0);
786 
787  permutedFFT(bStore, &*bStoreVecFreqPtr_);
788  permutedFFT(bstoreLeadCurrQ, &*bStoreLeadCurrQVecFreqPtr_);
789 
790  bLeadCurrentVecFreqPtr_->putScalar(0.0);
791  bLeadCurrentQVecFreqPtr_->putScalar(0.0);
792  permutedFFT(bNextLeadF, &*bLeadCurrentVecFreqPtr_);
794 
795 // Xyce::dout() << "HB Store Vector FD" << std::endl;
796 // bStoreVecFreqPtr_->printPetraObject(std::cout);
797 // bStoreLeadCurrQVecFreqPtr_->printPetraObject(std::cout);
798 // bLeadCurrentVecFreqPtr_->printPetraObject(Xyce::dout());
799 
800  int blockCount = bXf.blockCount();
801  int blockSize = bXf.block(0).globalLength();
802 // double omega = 2.0 * M_PI/ period_;
803 
804  int size_ = freqs_.size();
805  int posFreq = (size_-1)/2;
806  double omega = 2.0 * M_PI * freqs_[posFreq];
807 
808  for( int i = 0; i < blockCount; ++i )
809  {
810  // Create work vectors from the current frequency block vector
811  // NOTE: This needs to be done for each block to make sure that the
812  // map is the same as the bF block.
813  Linear::Vector QVec(bQ->block(i));
814  Linear::Vector freqVec = bQ->block(i);
815 
816  Linear::Vector dQdxdVpVec(bdQdxdVp->block(i));
817  Linear::Vector freqVec1 = bdQdxdVp->block(i);
818 
819  // Only one processor owns each block of the frequency-domain vector
820  if (freqVec.localLength() > 0)
821  {
822  omega = 2.0 * M_PI * freqs_[posFreq];
823 
824  QVec[0] = -freqVec[1]*omega;
825  QVec[1] = freqVec[0]*omega;
826 
827  dQdxdVpVec[0] = -freqVec1[1]*omega;
828  dQdxdVpVec[1] = freqVec1[0]*omega;
829 
830  for (int j=1; j < (blockSize/2+1)/2; ++j)
831  {
832 
833  omega = 2.0 * M_PI * freqs_[posFreq + j];
834 
835  QVec[2*j] = -freqVec[2*j+1]*omega;
836  QVec[2*(blockSize/2-j)] = -freqVec[2*j+1]*omega;
837 
838  QVec[2*j+1] = freqVec[2*j]*omega;
839  QVec[2*(blockSize/2-j)+1] = -freqVec[2*j]*omega;
840 
841  dQdxdVpVec[2*j] = -freqVec1[2*j+1]*omega;
842  dQdxdVpVec[2*(blockSize/2-j)] = -freqVec1[2*j+1]*omega;
843 
844  dQdxdVpVec[2*j+1] = freqVec1[2*j]*omega;
845  dQdxdVpVec[2*(blockSize/2-j)+1] = -freqVec1[2*j]*omega;
846 
847  }
848  }
849 
850  bF->block(i).update(1.0, QVec , 1.0);
851 
852  bdFdxdVp->block(i).update(1.0, dQdxdVpVec, 1.0);
853  }
854 
855  blockCount = bStoreVecFreqPtr_->blockCount();
856  blockSize = bStoreVecFreqPtr_->blockSize();
857 
858  for( int i = 0; i < blockCount; ++i )
859  {
860  // Create work vectors from the current frequency block vector
861  // NOTE: This needs to be done for each block to make sure that the
862  // map is the same as the bF block.
863  Linear::Vector stoLeadCurrdQdtVec(bStoreLeadCurrQVecFreqPtr_->block(i));
864  Linear::Vector stoLeadCurrQVec = bStoreLeadCurrQVecFreqPtr_->block(i);
865 
866  // Only one processor owns each block of the frequency-domain vector
867  if (stoLeadCurrQVec.localLength() > 0)
868  {
869 
870  omega = 2.0 * M_PI * freqs_[posFreq];
871  stoLeadCurrdQdtVec[0] = stoLeadCurrQVec[1]*omega;
872  stoLeadCurrdQdtVec[1] = stoLeadCurrQVec[0]*omega;
873 
874  for (int j=1; j < (blockSize/2+1)/2; ++j)
875  {
876  omega = 2.0 * M_PI * freqs_[posFreq+ j];
877 
878  stoLeadCurrdQdtVec[2*j] = -stoLeadCurrQVec[2*j+1]*omega;
879  stoLeadCurrdQdtVec[2*(blockSize/2-j)] = -stoLeadCurrQVec[2*j+1]*omega;
880 
881  stoLeadCurrdQdtVec[2*j+1] = stoLeadCurrQVec[2*j]*omega;
882  stoLeadCurrdQdtVec[2*(blockSize/2-j)+1] = -stoLeadCurrQVec[2*j]*omega;
883  }
884  }
885 
886  bStoreVecFreqPtr_->block(i).update(1.0, stoLeadCurrdQdtVec, 1.0);
887 // Xyce::dout() << "HB Store Vector dqdt + f(x) FD" << std::endl;
888 // bStoreVecFreqPtr_->printPetraObject(std:: cout);
889  }
890 
891  //
892  // take care of the lead current vector
893  blockCount = bLeadCurrentVecFreqPtr_->blockCount();
894  blockSize = bLeadCurrentVecFreqPtr_->blockSize();
895 
896  for( int i = 0; i < blockCount; ++i )
897  {
898  // Create work vectors from the current frequency block vector
899  // NOTE: This needs to be done for each block to make sure that the
900  // map is the same as the bF block.
901  Linear::Vector leadCurrdQdtVec(bLeadCurrentQVecFreqPtr_->block(i));
902  Linear::Vector leadCurrQVec = bLeadCurrentQVecFreqPtr_->block(i);
903 
904  // Only one processor owns each block of the frequency-domain vector
905  if (leadCurrQVec.localLength() > 0)
906  {
907 
908  omega = 2.0 * M_PI * freqs_[posFreq];
909  leadCurrdQdtVec[0] = leadCurrQVec[1]*omega;
910  leadCurrdQdtVec[1] = leadCurrQVec[0]*omega;
911 
912  for (int j=1; j < (blockSize/2+1)/2; ++j)
913  {
914  omega = 2.0 * M_PI * freqs_[posFreq+ j];
915 
916  leadCurrdQdtVec[2*j] = -leadCurrQVec[2*j+1]*omega;
917  leadCurrdQdtVec[2*(blockSize/2-j)] = -leadCurrQVec[2*j+1]*omega;
918 
919  leadCurrdQdtVec[2*j+1] = leadCurrQVec[2*j]*omega;
920  leadCurrdQdtVec[2*(blockSize/2-j)+1] = -leadCurrQVec[2*j]*omega;
921  }
922  }
923 
924  bLeadCurrentVecFreqPtr_->block(i).update(1.0, leadCurrdQdtVec, 1.0);
925  }
926 
927 
928 // permutedIFT(*bStoreVecFreqPtr_, &bStore);
929  if (DEBUG_HB)
930  {
931  Xyce::dout() << "HB X Vector" << std::endl;
932  bX.printPetraObject(std::cout);
933  // Xyce::dout() << "HB S Vector" << std::endl;
934  // bS.printPetraObject(Xyce::dout());
935  // Xyce::dout() << "HB dSdt Vector" << std::endl;
936  // bdSdt.printPetraObject(Xyce::dout());
937  Xyce::dout() << "HB Store Vector" << std::endl;
938  bStore.printPetraObject(std::cout);
939  Xyce::dout() << "HB Q Vector" << std::endl;
940  bQ->printPetraObject(std::cout);
941  Xyce::dout() << "HB F Vector" << std::endl;
942  bF->printPetraObject(std::cout);
943  Xyce::dout() << "HB bdFdxdVp Vector" << std::endl;
944  bdFdxdVp->printPetraObject(std::cout);
945  Xyce::dout() << "HB bdQdxdVp Vector" << std::endl;
946  bdQdxdVp->printPetraObject(std::cout);
947 
948 #ifndef Xyce_FLEXIBLE_DAE_LOADS
949  Xyce::dout() << "HB bmdQdx_" << std::endl;
950  bmdQdxPtr_->printPetraObject(std::cout);
951  Xyce::dout() << "HB bmdFdx_" << std::endl;
952  bmdFdxPtr_->printPetraObject(std::cout);
953 #endif // Xyce_FLEXIBLE_DAE_LOADS
954  Xyce::dout() << Xyce::section_divider << std::endl;
955  }
956 
957  return true;
958 }
959 
960 //-----------------------------------------------------------------------------
961 // Function : HBLoader::loadDeviceErrorWeightMask
962 // Purpose :
963 // Special Notes :
964 // Scope : public
965 // Creator : Eric Keiter
966 // Creation Date : 11/1/2014
967 //-----------------------------------------------------------------------------
968 bool
970  Linear::Vector * deviceMask) const
971 {
972 #if 0
973  Linear::BlockVector & bDevMask = *dynamic_cast<Linear::BlockVector*>(deviceMask);
974 
975 #if 0
976  Xyce::dout() << "HBLoader::loadDeviceErrorWeightMask. Original (nonblock) deviceMask.size = "
977  << appVecPtr_->globalLength() <<std::endl;
978 #endif
979 
980  appVecPtr_->putScalar(1.0);
981  bool returnValue = deviceManager_.loadErrorWeightMask(&*appVecPtr_);
982 
983  int blockCount = bDevMask.blockCount();
984  int blockSize = bDevMask.blockSize();
985 
986 #if 0
987  Xyce::dout() << "bDevMask.blockCount = "<< blockCount <<std::endl;
988  Xyce::dout() << "bDevMask.blockSize = "<< blockSize <<std::endl;
989  appVecPtr_->printPetraObject(Xyce::dout());
990 #endif
991 
992  //Teuchos::RCP<N_PDS_ParMap> baseMap = Teuchos::rcp_const_cast<N_PDS_ParMap>( hbBuilderPtr_->getBaseStoreMap() );
993 
994  for( int i = 0; i < blockCount; ++i )
995  {
996  // See if this variable is owned by the local processor.
997  // If so, this processor owns the entire j-th block of the vector
998  //int lid = baseMap->globalToLocalIndex( i );
999 
1000  Linear::Vector& localVecRef = bDevMask.block(i);
1001 
1002  for (int j=0;j<blockSize;++j)
1003  {
1004  localVecRef[j] = (*appVecPtr_)[i];
1005  }
1006  }
1007 
1008 #if 0
1009  bDevMask.printPetraObject(Xyce::dout());
1010 #endif
1011 
1012  return returnValue;
1013 #else
1014  return true;
1015 #endif
1016 }
1017 
1018 //-----------------------------------------------------------------------------
1019 // Function : HBLoader::permutedFFT
1020 // Purpose :
1021 // Special Notes :
1022 // Scope : public
1023 // Creator : Ting Mei
1024 // Creation Date : 09/05/08
1025 //---------------------------------------------------------------------------
1026 void HBLoader::permutedFFT(const Linear::BlockVector & xt, Linear::BlockVector * xf)
1027 {
1028  // Call the function to compute the permuted FFT from the block system helper functions.
1029  computePermutedDFT( *dftInterface_, xt, xf );
1030 }
1031 
1032 //-----------------------------------------------------------------------------
1033 // Function : HBLoader::permutedIFT
1034 // Purpose :
1035 // Special Notes :
1036 // Scope : public
1037 // Creator : Ting Mei
1038 // Creation Date : 09/05/08
1039 //---------------------------------------------------------------------------
1040 void HBLoader::permutedIFT(const Linear::BlockVector & xf, Linear::BlockVector * xt)
1041 {
1042  // Call the function to compute the permuted IFT from the block system helper functions.
1043  computePermutedIFT( *dftInterface_, xf, xt );
1044 }
1045 
1046 } // namespace Loader
1047 } // namespace Xyce
HBLoader(const Teuchos::RCP< const N_MPDE_Discretization > discPtr, Device::DeviceMgr &device_manager, Linear::Builder &builder)
Teuchos::RCP< Linear::Vector > appLastStaVecPtr_
Teuchos::RCP< Linear::Vector > appCurrStoVecPtr_
Teuchos::RCP< Linear::BlockMatrix > bmdFdxPtr_
Linear::Builder & builder_
std::vector< double > periodicTimes_
bool loadErrorWeightMask(Linear::Vector *deviceMaskPtr)
Pure virtual class to augment a linear system.
void setFastTimes(const std::vector< double > &times)
bool updateState(Linear::Vector *nextSolVectorPtr, Linear::Vector *currSolVectorPtr, Linear::Vector *lastSolVectorPtr, Linear::Vector *nextStaVectorPtr, Linear::Vector *currStaVectorPtr, Linear::Vector *lastStaVectorPtr, Linear::Vector *nextStoVectorPtr, Linear::Vector *currStoVectorPtr, Linear::Vector *lastStoVectorPtr)
void setFastTime(double timeVal)
Device::DeviceMgr & deviceManager_
bool loadDAEMatrices(Linear::Vector *X, Linear::Vector *S, Linear::Vector *dSdt, Linear::Vector *Store, Linear::Matrix *dQdx, Linear::Matrix *dFdx)
Teuchos::RCP< Loader > appLoaderPtr_
Actually a CktLoader.
Teuchos::RCP< const N_MPDE_Discretization > fastTimeDiscPtr_
Teuchos::RCP< Linear::Matrix > appdQdxPtr_
Teuchos::RCP< Linear::Vector > appStoLeadCurrQVecPtr_
std::vector< Teuchos::RCP< Linear::Matrix > > vecAppdFdxPtr_
void permutedFFT(const Linear::BlockVector &xt, Linear::BlockVector *xf)
Teuchos::RCP< Linear::Vector > appCurrStaVecPtr_
Teuchos::RCP< Linear::HBBuilder > hbBuilderPtr_
Device::DeviceMgr & deviceManager_
Definition: N_ANP_HB.C:2068
Teuchos::RCP< Linear::Vector > appCurrLeadFVecPtr_
Teuchos::RCP< Linear::BlockVector > bXtPtr_
Teuchos::RCP< Linear::BlockVector > bLeadCurrentVecFreqPtr_
Teuchos::RCP< Linear::BlockVector > bStoreVecFreqPtr_
Teuchos::RCP< Linear::BlockVector > bVtPtr_
Teuchos::RCP< Linear::Vector > appNextJunctionVVecPtr_
Teuchos::RCP< Linear::BlockVector > bStoreLeadCurrQVecFreqPtr_
Teuchos::RCP< Linear::BlockMatrix > bmdQdxPtr_
void setHBFreqs(const std::vector< double > &freqs)
void registerHBBuilder(Teuchos::RCP< Linear::HBBuilder > hbBuilderPtr)
Teuchos::RCP< N_UTL_DFTInterfaceDecl< std::vector< double > > > dftInterface_
Teuchos::RCP< Linear::Vector > appLeadQVecPtr_
Teuchos::RCP< Linear::Vector > appLastLeadFVecPtr_
bool loadTimeDepDAEMatrices(Linear::Vector *X, Linear::Vector *S, Linear::Vector *dSdt, Linear::Vector *Store, Linear::Matrix *dQdx, Linear::Matrix *dFdx)
Teuchos::RCP< Linear::Matrix > appdFdxPtr_
bool applyDAEMatrices(Linear::Vector *X, Linear::Vector *S, Linear::Vector *dSdt, Linear::Vector *Store, const Linear::Vector &V, Linear::Vector *dQdxV, Linear::Vector *dFdxV)
std::vector< double > times_
Linear::Builder & builder_
Definition: N_ANP_HB.C:2069
#define M_PI
Teuchos::RCP< Linear::Vector > appLastStoVecPtr_
Teuchos::RCP< Linear::BlockVector > bLeadCurrentQVecFreqPtr_
Teuchos::RCP< Linear::Vector > appNextLeadFVecPtr_
Teuchos::RCP< Linear::Vector > appLastJunctionVVecPtr_
std::vector< double > freqs_
Teuchos::RCP< Linear::Vector > appVecPtr_
void permutedIFT(const Linear::BlockVector &xf, Linear::BlockVector *xt)
bool loadDAEVectors(Linear::Vector *X, Linear::Vector *currX, Linear::Vector *lastX, Linear::Vector *S, Linear::Vector *currS, Linear::Vector *lastS, Linear::Vector *dSdt, Linear::Vector *Store, Linear::Vector *currStore, Linear::Vector *lastStore, Linear::Vector *storeLeadCurrQ, Linear::Vector *nextLeadFVectorPtr, Linear::Vector *currLeadFVectorPtr, Linear::Vector *lastLeadFVectorPtr, Linear::Vector *nextLeadQVectorPtr, Linear::Vector *nextJunctionVVectorPtr, Linear::Vector *currJunctionVVectorPtr, Linear::Vector *lastJunctionVVectorPtr, Linear::Vector *Q, Linear::Vector *F, Linear::Vector *B, Linear::Vector *dFdxdVp, Linear::Vector *dQdxdVp)
Teuchos::RCP< Linear::Vector > appCurrJunctionVVecPtr_
std::vector< Teuchos::RCP< Linear::Matrix > > vecAppdQdxPtr_
Teuchos::RCP< Linear::Vector > appNextStoVecPtr_
Teuchos::RCP< Linear::Vector > appNextStaVecPtr_
bool loadDeviceErrorWeightMask(Linear::Vector *deviceMask) const