Xyce  6.1
N_NLS_NOX_AugmentLinSys_GStepping.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_GStepping.C,v $
27 //
28 // Purpose : Algorithm for augmenting the Jacobian for pseudo
29 // transient solves using vnode conductance.
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.16 $
41 //
42 // Revision Date : $Date: 2015/04/08 19:18:29 $
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 : GStepping::GStepping
65 // Purpose : constructor
66 // Special Notes :
67 // Scope : public
68 // Creator : Roger Pawlowski, SNL 9233
69 // Creation Date :
70 //-----------------------------------------------------------------------------
72  const std::vector<int>& vnodeGIDVec,
73  Linear::Vector* cloneVector,
74  double scaledEndValue,
75  double resCond) :
76  node_list_type_(NLT_VoltageNodes),
77  vnodeGIDVec_(vnodeGIDVec),
78  tmp_vector_ptr_(0),
79  scaled_end_value_(scaledEndValue),
80  residualConductance_(resCond)
81 {
82  tmp_vector_ptr_ = new Linear::Vector(*cloneVector);
83 }
84 
85 //-----------------------------------------------------------------------------
86 // Function : GStepping::GStepping
87 // Purpose : constructor
88 // Special Notes :
89 // Scope : public
90 // Creator : Roger Pawlowski, SNL 9233
91 // Creation Date :
92 //-----------------------------------------------------------------------------
94 GStepping(const Teuchos::RCP<Epetra_MapColoring>& color_map,
95  Linear::Vector* cloneVector,
96  double scaledEndValue,
97  double resCond) :
98  node_list_type_(NLT_AllVoltageUnknowns),
99  scaled_end_value_(scaledEndValue),
100  residualConductance_(resCond)
101 {
102  color_map_ = color_map;
103  tmp_vector_ptr_ = new Linear::Vector(*cloneVector);
104 }
105 
106 //-----------------------------------------------------------------------------
107 // Function : GStepping::~GStepping
108 // Purpose : destructor
109 // Special Notes :
110 // Scope : public
111 // Creator : Roger Pawlowski, SNL 9233
112 // Creation Date :
113 //-----------------------------------------------------------------------------
115 {
116  delete tmp_vector_ptr_;
117 }
118 
119 
120 //-----------------------------------------------------------------------------
121 // Function : GStepping::setProgressVariable
122 // Purpose :
123 // Special Notes :
124 // Scope : public
125 // Creator : Roger Pawlowski, SNL 9233
126 // Creation Date :
127 //-----------------------------------------------------------------------------
128 void GStepping::setProgressVariable(double conductance)
129 {
130  // Exponential Continuation (con param goes from +4 -> -log10(endValue))
131  conductance_ = pow(10.0, conductance) - pow(10.0, scaled_end_value_) + residualConductance_;
132 }
133 
134 //-----------------------------------------------------------------------------
135 // Function : GStepping::augmentResidual
136 // Purpose :
137 // Special Notes :
138 // Scope : public
139 // Creator : Roger Pawlowski, SNL 9233
140 // Creation Date :
141 //-----------------------------------------------------------------------------
142 void GStepping::augmentResidual(const Linear::Vector * solution,
143  Linear::Vector * residualVector)
144 {
146  {
147  std::vector<int>::const_iterator i = vnodeGIDVec_.begin();
148  std::vector<int>::const_iterator stop = vnodeGIDVec_.end();
149  for ( ; i < stop; ++i)
150  {
151  double value = conductance_ *
152  (const_cast<Linear::Vector*>(solution))->getElementByGlobalIndex(*i);
153 
154  residualVector->sumElementByGlobalIndex(*i, value);
155  }
156  }
157  else
158  {
159  for (std::size_t i = 0; i < tmp_vector_ptr_->localLength(); ++i)
160  {
161  if ( (*color_map_)[i] == 0)
162  {
163  (*residualVector)[i] += conductance_ * (const_cast<Linear::Vector&>(*solution))[i];
164  }
165  }
166  }
167 }
168 
169 
170 //-----------------------------------------------------------------------------
171 // Function : GStepping::augmentJacobian
172 // Purpose :
173 // Special Notes :
174 // Scope : public
175 // Creator : Roger Pawlowski, SNL 9233
176 // Creation Date :
177 //-----------------------------------------------------------------------------
178 void GStepping::augmentJacobian(Linear::Matrix * jacobian)
179 {
180  jacobian->getDiagonal(*tmp_vector_ptr_);
181 
183  {
184  std::vector<int>::const_iterator i = vnodeGIDVec_.begin();
185  std::vector<int>::const_iterator stop = vnodeGIDVec_.end();
186  for ( ; i < stop; ++i)
187  {
188  tmp_vector_ptr_->sumElementByGlobalIndex(*i, conductance_);
189  }
190  }
191  else
192  {
193  for (std::size_t i = 0; i < tmp_vector_ptr_->localLength(); ++i)
194  {
195  if ( (*color_map_)[i] == 0)
196  {
197  (*tmp_vector_ptr_)[i] += conductance_;
198  }
199  }
200  }
201 
202  jacobian->replaceDiagonal(*tmp_vector_ptr_);
203 }
204 
205 }}}
void setProgressVariable(double time_step_size)
Set the progress variable (time step size for pseudo transient).
Pure virtual class to augment a linear system.
double scaled_end_value_
low end of the exponential term.
const T & value(const ParameterBase &entity, const Descriptor &descriptor)
Returns the value of the parameter for the entity.
Definition: N_DEV_Pars.h:1224
GStepping(const std::vector< int > &vnodeGIDVec, Xyce::Linear::Vector *cloneVector, double endValue, double residCond=0)
Ctor for the voltage nodes as a GID list.
void augmentResidual(const Xyce::Linear::Vector *solution, Xyce::Linear::Vector *residual_vector)
Augments the Residual.
void augmentJacobian(Xyce::Linear::Matrix *jacobian)
Augments the Jacobian.
double residualConductance_
residual value of the conductance. Should almost always be zero
const std::vector< int > vnodeGIDVec_
List of voltage node GIDs.
NodeListType node_list_type_
Type of list we are using.
Teuchos::RCP< Epetra_MapColoring > color_map_
Color 0 are the voltage unknowns.
Xyce::Linear::Vector * tmp_vector_ptr_
Temporary vector used to store diagonal.