libMesh::Problem_Interface Class Reference
Inheritance diagram for libMesh::Problem_Interface:

Public Member Functions

 Problem_Interface (NoxNonlinearSolver< Number > *solver)
 
 ~Problem_Interface ()
 
bool computeF (const Epetra_Vector &x, Epetra_Vector &FVec, NOX::Epetra::Interface::Required::FillType fillType)
 Compute and return F. More...
 
bool computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac)
 Compute an explicit Jacobian. More...
 
bool computePrecMatrix (const Epetra_Vector &x, Epetra_RowMatrix &M)
 Compute the Epetra_RowMatrix M, that will be used by the Aztec preconditioner instead of the Jacobian. This is used when there is no explicit Jacobian present (i.e. Matrix-Free Newton-Krylov). This MUST BE an Epetra_RowMatrix since the Aztec preconditioners need to know the sparsity pattern of the matrix. Returns true if computation was successful. More...
 
bool computePreconditioner (const Epetra_Vector &x, Epetra_Operator &Prec, Teuchos::ParameterList *p)
 Computes a user supplied preconditioner based on input vector x. Returns true if computation was successful. More...
 

Public Attributes

NoxNonlinearSolver< Number > * _solver
 

Detailed Description

Definition at line 53 of file trilinos_nox_nonlinear_solver.C.

Constructor & Destructor Documentation

libMesh::Problem_Interface::Problem_Interface ( NoxNonlinearSolver< Number > *  solver)
explicit

Definition at line 81 of file trilinos_nox_nonlinear_solver.C.

81  :
82  _solver(solver)
83 { }
libMesh::Problem_Interface::~Problem_Interface ( )

Definition at line 85 of file trilinos_nox_nonlinear_solver.C.

86 { }

Member Function Documentation

bool libMesh::Problem_Interface::computeF ( const Epetra_Vector &  x,
Epetra_Vector &  FVec,
NOX::Epetra::Interface::Required::FillType  fillType 
)

Compute and return F.

Definition at line 88 of file trilinos_nox_nonlinear_solver.C.

References _solver, libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::DofMap::enforce_constraints_exactly(), libMesh::err, libMesh::System::get_dof_map(), libMesh::NonlinearSolver< T >::matvec, libMesh::NonlinearImplicitSystem::ComputeResidual::residual(), libMesh::NonlinearSolver< T >::residual, libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian::residual_and_jacobian(), libMesh::NonlinearSolver< T >::residual_and_jacobian_object, libMesh::NonlinearSolver< T >::residual_object, libMesh::ExplicitSystem::rhs, libMesh::System::solution, libMesh::START_LOG(), libMesh::STOP_LOG(), libMesh::EpetraVector< T >::swap(), libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

90 {
91  START_LOG("residual()", "TrilinosNoxNonlinearSolver");
92 
93  NonlinearImplicitSystem &sys = _solver->system();
94 
95  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm()), R(r, sys.comm());
96  EpetraVector<Number>& X_sys = *libmesh_cast_ptr<EpetraVector<Number>*>(sys.solution.get());
97  EpetraVector<Number>& R_sys = *libmesh_cast_ptr<EpetraVector<Number>*>(sys.rhs);
98 
99  // Use the systems update() to get a good local version of the parallel solution
100  X_global.swap(X_sys);
101  R.swap(R_sys);
102 
104  sys.update();
105 
106  // Swap back
107  X_global.swap(X_sys);
108  R.swap(R_sys);
109 
110  R.zero();
111 
112  //-----------------------------------------------------------------------------
113  // if the user has provided both function pointers and objects only the pointer
114  // will be used, so catch that as an error
115 
117  {
118  libMesh::err << "ERROR: cannot specifiy both a function and object to compute the Residual!" << std::endl;
119  libmesh_error();
120  }
121 
123  {
124  libMesh::err << "ERROR: cannot specifiy both a function and object to compute the combined Residual & Jacobian!" << std::endl;
125  libmesh_error();
126  }
127  //-----------------------------------------------------------------------------
128 
129  if (_solver->residual != NULL) _solver->residual (*sys.current_local_solution.get(), R, sys);
130  else if (_solver->residual_object != NULL) _solver->residual_object->residual (*sys.current_local_solution.get(), R, sys);
131  else if (_solver->matvec != NULL) _solver->matvec (*sys.current_local_solution.get(), &R, NULL, sys);
132  else if (_solver->residual_and_jacobian_object != NULL) _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), &R, NULL, sys);
133  else return false;
134 
135  R.close();
136  X_global.close();
137 
138  STOP_LOG("residual()", "TrilinosNoxNonlinearSolver");
139 
140  return true;
141 }
bool libMesh::Problem_Interface::computeJacobian ( const Epetra_Vector &  x,
Epetra_Operator &  Jac 
)

Compute an explicit Jacobian.

Definition at line 143 of file trilinos_nox_nonlinear_solver.C.

References _solver, libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::DofMap::enforce_constraints_exactly(), libMesh::err, libMesh::System::get_dof_map(), libMesh::NonlinearImplicitSystem::ComputeJacobian::jacobian(), libMesh::NonlinearSolver< T >::jacobian, libMesh::NonlinearSolver< T >::jacobian_object, libMesh::NonlinearSolver< T >::matvec, libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian::residual_and_jacobian(), libMesh::NonlinearSolver< T >::residual_and_jacobian_object, libMesh::System::solution, libMesh::START_LOG(), libMesh::STOP_LOG(), libMesh::EpetraVector< T >::swap(), libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

145 {
146  START_LOG("jacobian()", "TrilinosNoxNonlinearSolver");
147 
148  NonlinearImplicitSystem &sys = _solver->system();
149 
150  EpetraMatrix<Number> Jac(&dynamic_cast<Epetra_FECrsMatrix &>(jac), sys.comm());
151  EpetraVector<Number>& X_sys = *libmesh_cast_ptr<EpetraVector<Number>*>(sys.solution.get());
152  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm());
153 
154  // Set the dof maps
155  Jac.attach_dof_map(sys.get_dof_map());
156 
157  // Use the systems update() to get a good local version of the parallel solution
158  X_global.swap(X_sys);
159 
161  sys.update();
162 
163  X_global.swap(X_sys);
164 
165  //-----------------------------------------------------------------------------
166  // if the user has provided both function pointers and objects only the pointer
167  // will be used, so catch that as an error
169  {
170  libMesh::err << "ERROR: cannot specify both a function and object to compute the Jacobian!" << std::endl;
171  libmesh_error();
172  }
173 
175  {
176  libMesh::err << "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!" << std::endl;
177  libmesh_error();
178  }
179  //-----------------------------------------------------------------------------
180 
181  if (_solver->jacobian != NULL) _solver->jacobian (*sys.current_local_solution.get(), Jac, sys);
182  else if (_solver->jacobian_object != NULL) _solver->jacobian_object->jacobian (*sys.current_local_solution.get(), Jac, sys);
183  else if (_solver->matvec != NULL) _solver->matvec (*sys.current_local_solution.get(), NULL, &Jac, sys);
184  else if (_solver->residual_and_jacobian_object != NULL) _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), NULL, &Jac, sys);
185  else libmesh_error();
186 
187  Jac.close();
188  X_global.close();
189 
190  STOP_LOG("jacobian()", "TrilinosNoxNonlinearSolver");
191 
192  return true;
193 }
bool libMesh::Problem_Interface::computePrecMatrix ( const Epetra_Vector &  x,
Epetra_RowMatrix &  M 
)

Compute the Epetra_RowMatrix M, that will be used by the Aztec preconditioner instead of the Jacobian. This is used when there is no explicit Jacobian present (i.e. Matrix-Free Newton-Krylov). This MUST BE an Epetra_RowMatrix since the Aztec preconditioners need to know the sparsity pattern of the matrix. Returns true if computation was successful.

Definition at line 195 of file trilinos_nox_nonlinear_solver.C.

196 {
197 // libMesh::out << "ERROR: Problem_Interface::preconditionVector() - Use Explicit Jacobian only for this test problem!" << endl;
198  throw 1;
199 }
bool libMesh::Problem_Interface::computePreconditioner ( const Epetra_Vector &  x,
Epetra_Operator &  Prec,
Teuchos::ParameterList *  p 
)

Computes a user supplied preconditioner based on input vector x. Returns true if computation was successful.

Definition at line 201 of file trilinos_nox_nonlinear_solver.C.

References _solver, libMesh::ParallelObject::comm(), libMesh::TrilinosPreconditioner< T >::compute(), libMesh::System::current_local_solution, libMesh::DofMap::enforce_constraints_exactly(), libMesh::err, libMesh::System::get_dof_map(), libMesh::NonlinearImplicitSystem::ComputeJacobian::jacobian(), libMesh::NonlinearSolver< T >::jacobian, libMesh::NonlinearSolver< T >::jacobian_object, libMesh::TrilinosPreconditioner< T >::mat(), libMesh::NonlinearSolver< T >::matvec, libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian::residual_and_jacobian(), libMesh::NonlinearSolver< T >::residual_and_jacobian_object, libMesh::System::solution, libMesh::START_LOG(), libMesh::STOP_LOG(), libMesh::EpetraVector< T >::swap(), libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

204 {
205  START_LOG("preconditioner()", "TrilinosNoxNonlinearSolver");
206 
207  NonlinearImplicitSystem &sys = _solver->system();
208  TrilinosPreconditioner<Number> & tpc = dynamic_cast<TrilinosPreconditioner<Number> &>(prec);
209 
210  EpetraMatrix<Number> Jac(dynamic_cast<Epetra_FECrsMatrix *>(tpc.mat()),sys.comm());
211  EpetraVector<Number>& X_sys = *libmesh_cast_ptr<EpetraVector<Number>*>(sys.solution.get());
212  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm());
213 
214  // Set the dof maps
215  Jac.attach_dof_map(sys.get_dof_map());
216 
217  // Use the systems update() to get a good local version of the parallel solution
218  X_global.swap(X_sys);
219 
221  sys.update();
222 
223  X_global.swap(X_sys);
224 
225  //-----------------------------------------------------------------------------
226  // if the user has provided both function pointers and objects only the pointer
227  // will be used, so catch that as an error
229  {
230  libMesh::err << "ERROR: cannot specify both a function and object to compute the Jacobian!" << std::endl;
231  libmesh_error();
232  }
233 
235  {
236  libMesh::err << "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!" << std::endl;
237  libmesh_error();
238  }
239  //-----------------------------------------------------------------------------
240 
241  if (_solver->jacobian != NULL) _solver->jacobian (*sys.current_local_solution.get(), Jac, sys);
242  else if (_solver->jacobian_object != NULL) _solver->jacobian_object->jacobian (*sys.current_local_solution.get(), Jac, sys);
243  else if (_solver->matvec != NULL) _solver->matvec (*sys.current_local_solution.get(), NULL, &Jac, sys);
244  else if (_solver->residual_and_jacobian_object != NULL) _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), NULL, &Jac, sys);
245  else libmesh_error();
246 
247  Jac.close();
248  X_global.close();
249 
250  tpc.compute();
251 
252  STOP_LOG("preconditioner()", "TrilinosNoxNonlinearSolver");
253 
254  return true;
255 }

Member Data Documentation

NoxNonlinearSolver<Number>* libMesh::Problem_Interface::_solver

The documentation for this class was generated from the following file:

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:58 UTC

Hosted By:
SourceForge.net Logo