Poster of Linux kernelThe best gift for a Linux geek
Problem_Interface

Problem_Interface

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

Problem_Interface -  

SYNOPSIS


 

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.
bool computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac)
Compute an explicit Jacobian.
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.
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.  

Public Attributes


NoxNonlinearSolver< Number > * _solver
 

Detailed Description

Definition at line 48 of file trilinos_nox_nonlinear_solver.C.  

Constructor & Destructor Documentation

 

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

Definition at line 74 of file trilinos_nox_nonlinear_solver.C.

                                                                        :
  _solver(solver)
{ }
 

Problem_Interface::~Problem_Interface ()

Definition at line 78 of file trilinos_nox_nonlinear_solver.C.

{ }
 

Member Function Documentation

 

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

Compute and return F.

Definition at line 81 of file trilinos_nox_nonlinear_solver.C.

References _solver, EpetraVector< T >::close(), System::current_local_solution, AutoPtr< Tp >::get(), NonlinearSolver< T >::residual, System::solution, EpetraVector< T >::swap(), NonlinearSolver< T >::system(), System::update(), and EpetraVector< T >::zero().

{
  NonlinearImplicitSystem &sys = _solver->system();

  EpetraVector<Number> X_global(*const_cast<Epetra_Vector *>(&x)), R(r);
  EpetraVector<Number>& X_sys = *libmesh_cast_ptr<EpetraVector<Number>*>(sys.solution.get());

  // Use the systems update() to get a good local version of the parallel solution
  X_global.swap(X_sys);
  sys.update();
  X_global.swap(X_sys);
  
  R.zero();
  
  if( _solver->residual == NULL )
    return false;
  _solver->residual (*sys.current_local_solution.get(), R);
  R.close();
  X_global.close();
  return true;
}
 

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

Compute an explicit Jacobian.

Definition at line 106 of file trilinos_nox_nonlinear_solver.C.

{
  throw 1;
  
  // return problem.evaluate(NOX::Epetra::Interface::Required::Jac, &x, 0, 0);
}
 

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

{
//   cout << 'ERROR: Problem_Interface::preconditionVector() - Use Explicit Jacobian only for this test problem!' << endl;
   throw 1;
}
 

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

{
//   cout << 'ERROR: Problem_Interface::preconditionVector() - Use Explicit Jacobian only for this test problem!' << endl;
   throw 1;
}
 

Member Data Documentation

 

NoxNonlinearSolver<Number>* Problem_Interface::_solver

Definition at line 70 of file trilinos_nox_nonlinear_solver.C.

Referenced by computeF().

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Public Member Functions
Public Attributes
Detailed Description
Constructor & Destructor Documentation
Problem_Interface::Problem_Interface (NoxNonlinearSolver< Number > *solver)
Problem_Interface::~Problem_Interface ()
Member Function Documentation
bool Problem_Interface::computeF (const Epetra_Vector &x, Epetra_Vector &FVec, NOX::Epetra::Interface::Required::FillTypefillType)
bool Problem_Interface::computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac)
bool Problem_Interface::computePrecMatrix (const Epetra_Vector &x, Epetra_RowMatrix &M)
bool Problem_Interface::computePreconditioner (const Epetra_Vector &x, Epetra_Operator &Prec, Teuchos::ParameterList *p)
Member Data Documentation
NoxNonlinearSolver<Number>* Problem_Interface::_solver
Author

This document was created by man2html, using the manual pages.
Time: 21:53:05 GMT, April 16, 2011