libMesh
Public Member Functions | Public Attributes | List of all members
libMesh::Problem_Interface Class Reference
Inheritance diagram for libMesh::Problem_Interface:
[legend]

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. More...
 
bool computePreconditioner (const Epetra_Vector &x, Epetra_Operator &Prec, Teuchos::ParameterList *p)
 Computes a user supplied preconditioner based on input vector x. More...
 

Public Attributes

NoxNonlinearSolver< Number > * _solver
 

Detailed Description

Definition at line 45 of file trilinos_nox_nonlinear_solver.C.

Constructor & Destructor Documentation

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

Definition at line 84 of file trilinos_nox_nonlinear_solver.C.

Referenced by libMesh::NoxNonlinearSolver< T >::init().

84  :
85  _solver(solver)
86 {}
NoxNonlinearSolver< Number > * _solver
libMesh::Problem_Interface::~Problem_Interface ( )

Definition at line 90 of file trilinos_nox_nonlinear_solver.C.

91 {}

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 95 of file trilinos_nox_nonlinear_solver.C.

References _solver, libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_dof_map(), libmesh_nullptr, 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::sys, libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

98 {
99  LOG_SCOPE("residual()", "TrilinosNoxNonlinearSolver");
100 
101  NonlinearImplicitSystem & sys = _solver->system();
102 
103  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm()), R(r, sys.comm());
104  EpetraVector<Number> & X_sys = *cast_ptr<EpetraVector<Number> *>(sys.solution.get());
105  EpetraVector<Number> & R_sys = *cast_ptr<EpetraVector<Number> *>(sys.rhs);
106 
107  // Use the systems update() to get a good local version of the parallel solution
108  X_global.swap(X_sys);
109  R.swap(R_sys);
110 
112  sys.update();
113 
114  // Swap back
115  X_global.swap(X_sys);
116  R.swap(R_sys);
117 
118  R.zero();
119 
120  //-----------------------------------------------------------------------------
121  // if the user has provided both function pointers and objects only the pointer
122  // will be used, so catch that as an error
123 
125  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the Residual!");
126 
128  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
129 
131  _solver->residual(*sys.current_local_solution.get(), R, sys);
132 
134  _solver->residual_object->residual(*sys.current_local_solution.get(), R, sys);
135 
136  else if (_solver->matvec != libmesh_nullptr)
137  _solver->matvec(*sys.current_local_solution.get(), &R, libmesh_nullptr, sys);
138 
139  else if (_solver->residual_and_jacobian_object != libmesh_nullptr)
140  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), &R, libmesh_nullptr, sys);
141 
142  else
143  return false;
144 
145  R.close();
146  X_global.close();
147 
148  return true;
149 }
void(* residual)(const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)
Function that computes the residual R(X) of the nonlinear system at the input iterate X...
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
Object that computes either the residual or the Jacobian of the nonlinear system at the input itera...
virtual void residual(const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)=0
Residual function.
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
NoxNonlinearSolver< Number > * _solver
NonlinearImplicitSystem::ComputeResidual * residual_object
Object that computes the residual R(X) of the nonlinear system at the input iterate X...
PetscErrorCode Vec x
const DofMap & get_dof_map() const
Definition: system.h:2030
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=libmesh_nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Definition: dof_map.h:1816
const sys_type & system() const
virtual void residual_and_jacobian(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)=0
Residual & Jacobian function, calculated simultaneously.
void(* matvec)(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
Function that computes either the residual or the Jacobian of the nonlinear system at the input ite...
bool libMesh::Problem_Interface::computeJacobian ( const Epetra_Vector &  x,
Epetra_Operator &  Jac 
)

Compute an explicit Jacobian.

Definition at line 153 of file trilinos_nox_nonlinear_solver.C.

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

155 {
156  LOG_SCOPE("jacobian()", "TrilinosNoxNonlinearSolver");
157 
158  NonlinearImplicitSystem & sys = _solver->system();
159 
160  EpetraMatrix<Number> Jac(&dynamic_cast<Epetra_FECrsMatrix &>(jac), sys.comm());
161  EpetraVector<Number> & X_sys = *cast_ptr<EpetraVector<Number> *>(sys.solution.get());
162  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm());
163 
164  // Set the dof maps
165  Jac.attach_dof_map(sys.get_dof_map());
166 
167  // Use the systems update() to get a good local version of the parallel solution
168  X_global.swap(X_sys);
169 
171  sys.update();
172 
173  X_global.swap(X_sys);
174 
175  //-----------------------------------------------------------------------------
176  // if the user has provided both function pointers and objects only the pointer
177  // will be used, so catch that as an error
179  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the Jacobian!");
180 
182  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
183 
185  _solver->jacobian(*sys.current_local_solution.get(), Jac, sys);
186 
188  _solver->jacobian_object->jacobian(*sys.current_local_solution.get(), Jac, sys);
189 
190  else if (_solver->matvec != libmesh_nullptr)
191  _solver->matvec(*sys.current_local_solution.get(), libmesh_nullptr, &Jac, sys);
192 
194  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), libmesh_nullptr, &Jac, sys);
195 
196  else
197  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
198 
199  Jac.close();
200  X_global.close();
201 
202  return true;
203 }
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
Object that computes either the residual or the Jacobian of the nonlinear system at the input itera...
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
NoxNonlinearSolver< Number > * _solver
void(* jacobian)(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
PetscErrorCode Vec x
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual void jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)=0
Jacobian function.
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=libmesh_nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Definition: dof_map.h:1816
const sys_type & system() const
virtual void residual_and_jacobian(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)=0
Residual & Jacobian function, calculated simultaneously.
void(* matvec)(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
Function that computes either the residual or the Jacobian of the nonlinear system at the input ite...
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 207 of file trilinos_nox_nonlinear_solver.C.

208 {
209  // libMesh::out << "ERROR: Problem_Interface::preconditionVector() - Use Explicit Jacobian only for this test problem!" << endl;
210  throw 1;
211 }
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 215 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::System::get_dof_map(), libMesh::NonlinearImplicitSystem::ComputeJacobian::jacobian(), libMesh::NonlinearSolver< T >::jacobian, libMesh::NonlinearSolver< T >::jacobian_object, libmesh_nullptr, 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::sys, libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

218 {
219  LOG_SCOPE("preconditioner()", "TrilinosNoxNonlinearSolver");
220 
221  NonlinearImplicitSystem & sys = _solver->system();
222  TrilinosPreconditioner<Number> & tpc = dynamic_cast<TrilinosPreconditioner<Number> &>(prec);
223 
224  EpetraMatrix<Number> Jac(dynamic_cast<Epetra_FECrsMatrix *>(tpc.mat()),sys.comm());
225  EpetraVector<Number> & X_sys = *cast_ptr<EpetraVector<Number> *>(sys.solution.get());
226  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm());
227 
228  // Set the dof maps
229  Jac.attach_dof_map(sys.get_dof_map());
230 
231  // Use the systems update() to get a good local version of the parallel solution
232  X_global.swap(X_sys);
233 
235  sys.update();
236 
237  X_global.swap(X_sys);
238 
239  //-----------------------------------------------------------------------------
240  // if the user has provided both function pointers and objects only the pointer
241  // will be used, so catch that as an error
243  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the Jacobian!");
244 
246  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
247 
249  _solver->jacobian(*sys.current_local_solution.get(), Jac, sys);
250 
252  _solver->jacobian_object->jacobian(*sys.current_local_solution.get(), Jac, sys);
253 
254  else if (_solver->matvec != libmesh_nullptr)
255  _solver->matvec(*sys.current_local_solution.get(), libmesh_nullptr, &Jac, sys);
256 
258  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), libmesh_nullptr, &Jac, sys);
259 
260  else
261  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
262 
263  Jac.close();
264  X_global.close();
265 
266  tpc.compute();
267 
268  return true;
269 }
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
Object that computes either the residual or the Jacobian of the nonlinear system at the input itera...
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
NoxNonlinearSolver< Number > * _solver
void(* jacobian)(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
PetscErrorCode Vec x
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual void jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)=0
Jacobian function.
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=libmesh_nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Definition: dof_map.h:1816
const sys_type & system() const
virtual void residual_and_jacobian(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)=0
Residual & Jacobian function, calculated simultaneously.
void(* matvec)(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
Function that computes either the residual or the Jacobian of the nonlinear system at the input ite...

Member Data Documentation

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

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