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

◆ Problem_Interface()

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

Definition at line 84 of file trilinos_nox_nonlinear_solver.C.

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

◆ ~Problem_Interface()

libMesh::Problem_Interface::~Problem_Interface ( )
default

Member Function Documentation

◆ computeF()

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

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

97 {
98  LOG_SCOPE("residual()", "TrilinosNoxNonlinearSolver");
99 
100  NonlinearImplicitSystem & sys = _solver->system();
101 
102  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm()), R(r, sys.comm());
103  EpetraVector<Number> & X_sys = *cast_ptr<EpetraVector<Number> *>(sys.solution.get());
104  EpetraVector<Number> & R_sys = *cast_ptr<EpetraVector<Number> *>(sys.rhs);
105 
106  // Use the systems update() to get a good local version of the parallel solution
107  X_global.swap(X_sys);
108  R.swap(R_sys);
109 
111  sys.get_dof_map().enforce_constraints_exactly(sys);
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 
124  libmesh_error_msg_if(_solver->residual && _solver->residual_object, "ERROR: cannot specify both a function and object to compute the Residual!");
125 
126  libmesh_error_msg_if(_solver->matvec && _solver->residual_and_jacobian_object,
127  "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
128 
129  if (_solver->residual != nullptr)
130  _solver->residual(*sys.current_local_solution.get(), R, sys);
131 
132  else if (_solver->residual_object != nullptr)
133  _solver->residual_object->residual(*sys.current_local_solution.get(), R, sys);
134 
135  else if (_solver->matvec != nullptr)
136  _solver->matvec(*sys.current_local_solution.get(), &R, nullptr, sys);
137 
138  else if (_solver->residual_and_jacobian_object != nullptr)
139  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), &R, nullptr, sys);
140 
141  else
142  return false;
143 
144  R.close();
145  X_global.close();
146 
147  return true;
148 }
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.
NoxNonlinearSolver< Number > * _solver
template class LIBMESH_EXPORT EpetraVector< Number >
NonlinearImplicitSystem::ComputeResidual * residual_object
Object that computes the residual R(X) of the nonlinear system at the input iterate X...
const sys_type & system() const
bool _exact_constraint_enforcement
Whether we should enforce exact constraints globally during a solve.
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...

◆ computeJacobian()

bool libMesh::Problem_Interface::computeJacobian ( const Epetra_Vector &  x,
Epetra_Operator &  Jac 
)

Compute an explicit Jacobian.

Definition at line 152 of file trilinos_nox_nonlinear_solver.C.

References libMesh::NonlinearSolver< T >::_exact_constraint_enforcement, _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::NonlinearSolver< T >::matvec, libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian::residual_and_jacobian(), libMesh::NonlinearSolver< T >::residual_and_jacobian_object, libMesh::System::solution, libMesh::EpetraVector< T >::swap(), libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

154 {
155  LOG_SCOPE("jacobian()", "TrilinosNoxNonlinearSolver");
156 
157  NonlinearImplicitSystem & sys = _solver->system();
158 
159  EpetraMatrix<Number> J(&dynamic_cast<Epetra_FECrsMatrix &>(jac), sys.comm());
160  EpetraVector<Number> & X_sys = *cast_ptr<EpetraVector<Number> *>(sys.solution.get());
161  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x), sys.comm());
162 
163  // Set the dof maps
164  J.attach_dof_map(sys.get_dof_map());
165 
166  // Use the systems update() to get a good local version of the parallel solution
167  X_global.swap(X_sys);
168 
170  sys.get_dof_map().enforce_constraints_exactly(sys);
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
178  libmesh_error_msg_if(_solver->jacobian && _solver->jacobian_object,
179  "ERROR: cannot specify both a function and object to compute the Jacobian!");
180 
181  libmesh_error_msg_if(_solver->matvec && _solver->residual_and_jacobian_object,
182  "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
183 
184  if (_solver->jacobian != nullptr)
185  _solver->jacobian(*sys.current_local_solution.get(), J, sys);
186 
187  else if (_solver->jacobian_object != nullptr)
188  _solver->jacobian_object->jacobian(*sys.current_local_solution.get(), J, sys);
189 
190  else if (_solver->matvec != nullptr)
191  _solver->matvec(*sys.current_local_solution.get(), nullptr, &J, sys);
192 
193  else if (_solver->residual_and_jacobian_object != nullptr)
194  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), nullptr, &J, sys);
195 
196  else
197  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
198 
199  J.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...
NoxNonlinearSolver< Number > * _solver
template class LIBMESH_EXPORT EpetraVector< Number >
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...
const sys_type & system() const
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
bool _exact_constraint_enforcement
Whether we should enforce exact constraints globally during a solve.
virtual void jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)=0
Jacobian function.
template class LIBMESH_EXPORT EpetraMatrix< Number >
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...

◆ computePrecMatrix()

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 }

◆ computePreconditioner()

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::TrilinosPreconditioner< T >::mat(), libMesh::NonlinearSolver< T >::matvec, libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian::residual_and_jacobian(), libMesh::NonlinearSolver< T >::residual_and_jacobian_object, libMesh::System::solution, libMesh::EpetraVector< T >::swap(), libMesh::NonlinearSolver< T >::system(), and libMesh::System::update().

218 {
219  LOG_SCOPE("preconditioner()", "TrilinosNoxNonlinearSolver");
220 
221  NonlinearImplicitSystem & sys = _solver->system();
223 
224  EpetraMatrix<Number> J(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  J.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 
234  if (this->_exact_constraint_enforcement)
235  sys.get_dof_map().enforce_constraints_exactly(sys);
236  sys.update();
237 
238  X_global.swap(X_sys);
239 
240  //-----------------------------------------------------------------------------
241  // if the user has provided both function pointers and objects only the pointer
242  // will be used, so catch that as an error
243  libmesh_error_msg_if(_solver->jacobian && _solver->jacobian_object,
244  "ERROR: cannot specify both a function and object to compute the Jacobian!");
245 
246  libmesh_error_msg_if(_solver->matvec && _solver->residual_and_jacobian_object,
247  "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
248 
249  if (_solver->jacobian != nullptr)
250  _solver->jacobian(*sys.current_local_solution.get(), J, sys);
251 
252  else if (_solver->jacobian_object != nullptr)
253  _solver->jacobian_object->jacobian(*sys.current_local_solution.get(), J, sys);
254 
255  else if (_solver->matvec != nullptr)
256  _solver->matvec(*sys.current_local_solution.get(), nullptr, &J, sys);
257 
258  else if (_solver->residual_and_jacobian_object != nullptr)
259  _solver->residual_and_jacobian_object->residual_and_jacobian (*sys.current_local_solution.get(), nullptr, &J, sys);
260 
261  else
262  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
263 
264  J.close();
265  X_global.close();
266 
267  tpc.compute();
268 
269  return true;
270 }
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
Object that computes either the residual or the Jacobian of the nonlinear system at the input itera...
template class LIBMESH_EXPORT TrilinosPreconditioner< Number >
NoxNonlinearSolver< Number > * _solver
template class LIBMESH_EXPORT EpetraVector< Number >
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...
const sys_type & system() const
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
virtual void jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)=0
Jacobian function.
template class LIBMESH_EXPORT EpetraMatrix< Number >
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

◆ _solver

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

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