Xyce  6.1
N_NLS_NOX_SharedSystem.h
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_SharedSystem.h,v $
27 //
28 // Purpose : Interface to let multiple N_NLS::NOX::Group's
29 // share a single system of RHS Vector, Jacobian
30 // matrix, Newton vector, and gradient vector.
31 // Closely related to the NOX::Epetra::SharedSystem
32 // class.
33 //
34 // Special Notes :
35 //
36 // Creator : Tammy Kolda, NLS, 8950
37 //
38 // Creation Date : 01/31/02
39 //
40 // Revision Information:
41 // ---------------------
42 //
43 // Revision Number: $Revision: 1.29.2.1 $
44 //
45 // Revision Date : $Date: 2015/04/02 18:20:17 $
46 //
47 // Current Owner : $Author: tvrusso $
48 //-------------------------------------------------------------------------
49 
50 #ifndef Xyce_N_NLS_NOX_SharedSystem_h
51 #define Xyce_N_NLS_NOX_SharedSystem_h
52 
53 // ---------- Standard Includes ----------
54 
55 // ---------- Xyce Includes ----------
56 
57 // ---------- NOX Includes ----------
58 
59 #include <N_NLS_fwd.h>
60 
61 #include "N_NLS_NOX_Vector.h"
62 #include "N_NLS_NOX_Group.h"
63 
64 class Ifpack_IlukGraph;
65 class Ifpack_CrsRiluk;
66 
67 //-----------------------------------------------------------------------------
68 // Class : N_NLS::NOX::SharedSystem
69 //
70 // Purpose :
71 //
72 // Interface to let multiple N_NLS::NOX::Group's share the
73 // vectors and matrices in the Xyce nonlinear solver.
74 //
75 // Closely related conceptually to the
76 // NOX::Epetra::SharedJacobian class.
77 //
78 // Creator : Tammy Kolda, SNL, 8950
79 //
80 // Creation Date : 1/31/02
81 //-----------------------------------------------------------------------------
82 
83 namespace Xyce {
84 namespace Nonlinear {
85 namespace N_NLS_NOX {
86 class SharedSystem {
87 
88 public:
89 
90  SharedSystem(Xyce::Linear::Vector& x,
91  Xyce::Linear::Vector& f,
92  Xyce::Linear::Matrix& jacobian,
93  Xyce::Linear::Vector& newton,
94  Xyce::Linear::Vector& gradient,
95  Xyce::Linear::System& lasSys,
96  Interface& interface);
97 
98  ~SharedSystem();
99 
100 
101  void reset(Xyce::Linear::Vector& x,
102  Xyce::Linear::Vector& f,
103  Xyce::Linear::Matrix& jacobian,
104  Xyce::Linear::Vector& newton,
105  Xyce::Linear::Vector& gradient,
106  Xyce::Linear::System& lasSys,
107  Interface& interface);
108 
109  //---------------------------------------------------------------------------
110  // Function : isJacobianOwner
111  // Purpose : Verify that the group pointed to by grp is owner of the
112  // Jacobian matrix.
113  //---------------------------------------------------------------------------
114  inline bool isJacobianOwner(const Group* grp) const
115  {
116  return (grp == ownerOfJacobian_);
117  };
118 
119  //---------------------------------------------------------------------------
120  // Function : areStateVectors
121  // Purpose : To compute a Jacobian, the state vectors must be
122  // updated with respect to the solution in the group.
123  // However, the state vectors are updated ONLY during
124  // calls to compute the residual. This method checks
125  // to see if the state vectors still correspond to this
126  // group. Returns true if state vectors are correct.
127  //---------------------------------------------------------------------------
128  inline bool areStateVectors(const Group* grp) const
129  {
130  return (grp == ownerOfStateVectors_);
131  };
132 
133  bool computeF(const Vector& solution, Vector& F, const Group* grp);
134 
135  bool computeJacobian(Group* grp);
136 
137  bool computeNewton(const Vector& F, Vector& Newton,
138  Teuchos::ParameterList& params);
139 
140  bool computeGradient(const Vector& F, Vector& Gradient);
141 
142  bool applyJacobian(const Vector& input, Vector& result) const;
143 
144  bool applyJacobianTranspose(const Vector& input, Vector& result) const;
145 
147 
148  // Take ownership of const Jacobian.
149  const Xyce::Linear::Matrix& getJacobian() const;
150 
151  // Take ownership of Jacobian and get a reference to it.
152  Xyce::Linear::Matrix& getJacobian(const Group* grp);
153 
154  // Take ownership of the state vectors.
155  void getStateVectors(const Group* grp);
156 
157  // Get a pointer to the Xyce::Linear::System object
158  Xyce::Linear::System* getLasSystem();
159 
160  // Use for debugging (corresponding to the ones in N_NLS_NonLinearSolver).
161  void debugOutput1 (Xyce::Linear::Matrix & jacobian, Xyce::Linear::Vector & rhs);
162  void debugOutput3 (Xyce::Linear::Vector & dxVector, Xyce::Linear::Vector & xVector);
163 
164  // Preconditioning objects for the Group::applyRightPreconditioning method
165  bool computePreconditioner();
166  bool deletePreconditioner();
167  bool applyRightPreconditioning(bool useTranspose,
168  Teuchos::ParameterList& params,
169  const Vector& input,
170  Vector& result);
171 
172  // This is used to construct vectors in the group. We used to
173  // clone a time integrator vector but when the DC Op point fails,
174  // somewhere (I have no idea where) it decides to delete vectors
175  // before the Nonlinear solver can delete theirs. This causes
176  // a seg fault.
177  //
178  Vector* cloneSolutionVector() const;
179 
180  // Take ownership of const newton vector
181  const Vector & getNewtonVector() const;
182 
183  void printSoln(std::ostream &os) {xyceSolnPtr_->print(os);}
184  void printRes(std::ostream &os) {xyceFPtr_->print(os);}
185 
186 private:
187 
188  // Views of xyce objects used in the fill
189  Vector* xyceSolnPtr_; // Solution vector
190  Vector* xyceFPtr_; // Residual Vector
191  Xyce::Linear::Matrix* xyceJacobianPtr_; // Jacobian matrix
192  Vector* xyceNewtonPtr_; // Newton Vector
193  Vector* xyceGradientPtr_; // gradient Vector
194  Xyce::Linear::System* xyceLasSysPtr_; // LAS System
195  Interface* xyceInterfacePtr_; // Nonlinear Solver Interface
196 
197  // Flag for Matrix Free Loads tscoffe/tmei 07/29/08
199 
202 
203  // Ifpack preconditioning objects for applyRightPreconditioning method
204  mutable Ifpack_IlukGraph* ifpackGraphPtr_;
205  mutable Ifpack_CrsRiluk* ifpackPreconditionerPtr_;
206 
207 }; // class SharedSystem
208 }}} // namespace N_NLS_NOX
209 
210 #endif // Xyce_N_NLS_NOX_SharedSystem_h
211 
bool areStateVectors(const Group *grp) const
bool applyJacobian(const Vector &input, Vector &result) const
void reset(Xyce::Linear::Vector &x, Xyce::Linear::Vector &f, Xyce::Linear::Matrix &jacobian, Xyce::Linear::Vector &newton, Xyce::Linear::Vector &gradient, Xyce::Linear::System &lasSys, Interface &interface)
bool computeNewton(const Vector &F, Vector &Newton, Teuchos::ParameterList &params)
Pure virtual class to augment a linear system.
bool applyJacobianTranspose(const Vector &input, Vector &result) const
void debugOutput3(Xyce::Linear::Vector &dxVector, Xyce::Linear::Vector &xVector)
bool applyRightPreconditioning(bool useTranspose, Teuchos::ParameterList &params, const Vector &input, Vector &result)
SharedSystem(Xyce::Linear::Vector &x, Xyce::Linear::Vector &f, Xyce::Linear::Matrix &jacobian, Xyce::Linear::Vector &newton, Xyce::Linear::Vector &gradient, Xyce::Linear::System &lasSys, Interface &interface)
const Xyce::Linear::Matrix & getJacobian() const
bool computeGradient(const Vector &F, Vector &Gradient)
bool isJacobianOwner(const Group *grp) const
bool computeF(const Vector &solution, Vector &F, const Group *grp)
void debugOutput1(Xyce::Linear::Matrix &jacobian, Xyce::Linear::Vector &rhs)
void print(std::ostream &os) const