Xyce  6.1
N_NLS_NOX_AugmentLinSys_PseudoTransient.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_NLS_NOX_AugmentLinSys_PseudoTransient.C,v $
27 //
28 // Purpose : Algorithm for augmenting the Jacobian for pseudo
29 // transient solves.
30 //
31 // Special Notes :
32 //
33 // Creator : Roger Pawlowski, SNL 9233
34 //
35 // Creation Date : 03/07/06
36 //
37 // Revision Information:
38 // ---------------------
39 //
40 // Revision Number: $Revision: 1.11.2.1 $
41 //
42 // Revision Date : $Date: 2015/04/02 18:20:17 $
43 //
44 // Current Owner : $Author: tvrusso $
45 //-------------------------------------------------------------------------
46 
47 #include <Xyce_config.h>
48 
49 
50 // ---------- Standard Includes ----------
51 
52 // ---------- Xyce Includes ----------
53 
54 #include "N_LAS_Vector.h"
55 #include "N_LAS_Matrix.h"
56 #include "Epetra_MapColoring.h"
58 
59 namespace Xyce {
60 namespace Nonlinear {
61 namespace N_NLS_NOX {
62 
63 //-----------------------------------------------------------------------------
64 // Function : AugmentLinSysPseudoTransient::AugmentLinSysPseudoTransient
65 // Purpose : constructor
66 // Special Notes :
67 // Scope : public
68 // Creator : Roger Pawlowski, SNL 9233
69 // Creation Date :
70 //-----------------------------------------------------------------------------
72  (const Teuchos::RCP<Epetra_MapColoring>&
73  color_map,
74  Linear::Vector* cloneVector,
75  bool useVoltageScaleFactor,
76  double voltageScaleFactor)
77 {
78  use_voltage_scale_factor_ = useVoltageScaleFactor;
79  voltage_scale_factor_ = voltageScaleFactor;
80  color_map_ = color_map;
81  tmp_vector_ptr_ = new Linear::Vector(*cloneVector);
82 }
83 
84 //-----------------------------------------------------------------------------
85 // Function : AugmentLinSysPseudoTransient::~AugmentLinSysPseudoTransient
86 // Purpose : destructor
87 // Special Notes :
88 // Scope : public
89 // Creator : Roger Pawlowski, SNL 9233
90 // Creation Date :
91 //-----------------------------------------------------------------------------
93 {
94  delete tmp_vector_ptr_;
95 }
96 
97 
98 //-----------------------------------------------------------------------------
99 // Function : AugmentLinSysPseudoTransient::setProgressVariable
100 // Purpose :
101 // Special Notes :
102 // Scope : public
103 // Creator : Roger Pawlowski, SNL 9233
104 // Creation Date :
105 //-----------------------------------------------------------------------------
107  (double time_step_size)
108 {
109  time_step_size_ = time_step_size;
110 }
111 
112 //-----------------------------------------------------------------------------
113 // Function : AugmentLinSysPseudoTransient::augmentResidual
114 // Purpose :
115 // Special Notes : no-op for pseudo-transient.
116 // Scope : public
117 // Creator : Roger Pawlowski, SNL 9233
118 // Creation Date :
119 //-----------------------------------------------------------------------------
121  (const Linear::Vector * solution, Linear::Vector * residual_vector)
122 {
123  // Nothing to do for pseudo transient!!
124 }
125 
126 //-----------------------------------------------------------------------------
127 // Function : AugmentLinSysPseudoTransient::augmentJacobian
128 // Purpose :
129 // Special Notes :
130 // Scope : public
131 // Creator : Roger Pawlowski, SNL 9233
132 // Creation Date :
133 //-----------------------------------------------------------------------------
135  (Linear::Matrix * jacobian)
136 {
137  //cout << "Augmenting Jacobian for Pseudo Transient" << endl;
138  //cout << "Pseudo Trans Step Size = " << pseudoTransientTimeStep_ << endl;
139 
140  //color_map_->Print(cout);
141 
142  //jacobian->printPetraObject();
143 
144  //jacobian->scale(conParamValue);
145 
146  jacobian->getDiagonal(*tmp_vector_ptr_);
147 
148  //tmp_vector_ptr_->printPetraObject();
149 
150  double value = 1.0 / time_step_size_;
151 
152  //cout << "Pseudo Transient Time Step Size = " << time_step_size_ << endl;
153 
154  if (!use_voltage_scale_factor_)
155  {
156  tmp_vector_ptr_->addScalar(value);
157  }
158  else
159  {
160  for (std::size_t i = 0; i < tmp_vector_ptr_->localLength(); ++i)
161  {
162  if ( (*color_map_)[i] == 0)
163  {
164  (*tmp_vector_ptr_)[i] += value * voltage_scale_factor_;
165  }
166  else
167  {
168  (*tmp_vector_ptr_)[i] += value;
169  }
170  }
171  //RPP Might need to export local values for tmp_vector_ptr_ here
172  //for parallel.
173  }
174 
175  jacobian->replaceDiagonal(*tmp_vector_ptr_);
176 
177  //jacobian->printPetraObject();
178 
179  //cout << "Time step size = " << time_step_size_ << endl;
180 
181  //if (use_voltage_scale_factor_)
182  //cout << "Voltage scale factor = " << voltage_scale_factor_ << endl;
183 }
184 
185 }}}
Pure virtual class to augment a linear system.
const T & value(const ParameterBase &entity, const Descriptor &descriptor)
Returns the value of the parameter for the entity.
Definition: N_DEV_Pars.h:1224
void setProgressVariable(double time_step_size)
Set the progress variable (time step size for pseudo transient).
void augmentJacobian(Xyce::Linear::Matrix *jacobian)
Augments the Jacobian.
void augmentResidual(const Xyce::Linear::Vector *solution, Xyce::Linear::Vector *residual_vector)
Augments the Residual.
AugmentLinSysPseudoTransient(const Teuchos::RCP< Epetra_MapColoring > &color_map, Xyce::Linear::Vector *cloneVector, bool useVoltageScaleFactor=false, double voltageScaleFactor=1.0)
Ctor.