![]() |
|
#include <iterative_linear_solver.h>
Inheritance diagram for oomph::IterativeLinearSolver:Public Member Functions | |
| IterativeLinearSolver () | |
| IterativeLinearSolver (const IterativeLinearSolver &)=delete | |
| Broken copy constructor. More... | |
| void | operator= (const IterativeLinearSolver &)=delete |
| Broken assignment operator. More... | |
| virtual | ~IterativeLinearSolver () |
| Destructor (empty) More... | |
| Preconditioner *& | preconditioner_pt () |
| Access function to preconditioner. More... | |
| Preconditioner *const & | preconditioner_pt () const |
| Access function to preconditioner (const version) More... | |
| double & | tolerance () |
| Access to convergence tolerance. More... | |
| unsigned & | max_iter () |
| Access to max. number of iterations. More... | |
| virtual unsigned | iterations () const =0 |
| Number of iterations taken. More... | |
| void | enable_doc_convergence_history () |
| Enable documentation of the convergence history. More... | |
| void | disable_doc_convergence_history () |
| Disable documentation of the convergence history. More... | |
| void | open_convergence_history_file_stream (const std::string &file_name, const std::string &zone_title="") |
| void | close_convergence_history_file_stream () |
| Close convergence history output stream. More... | |
| double | jacobian_setup_time () const |
| double | linear_solver_solution_time () const |
| return the time taken to solve the linear system More... | |
| virtual double | preconditioner_setup_time () const |
| returns the the time taken to setup the preconditioner More... | |
| void | enable_setup_preconditioner_before_solve () |
| Setup the preconditioner before the solve. More... | |
| void | disable_setup_preconditioner_before_solve () |
| Don't set up the preconditioner before the solve. More... | |
| void | enable_error_after_max_iter () |
| Throw an error if we don't converge within max_iter. More... | |
| void | disable_error_after_max_iter () |
| Don't throw an error if we don't converge within max_iter (default). More... | |
| void | enable_iterative_solver_as_preconditioner () |
| void | disable_iterative_solver_as_preconditioner () |
Public Member Functions inherited from oomph::LinearSolver | |
| LinearSolver () | |
| Empty constructor, initialise the member data. More... | |
| LinearSolver (const LinearSolver &dummy)=delete | |
| Broken copy constructor. More... | |
| void | operator= (const LinearSolver &)=delete |
| Broken assignment operator. More... | |
| virtual | ~LinearSolver () |
| Empty virtual destructor. More... | |
| void | enable_doc_time () |
| Enable documentation of solve times. More... | |
| void | disable_doc_time () |
| Disable documentation of solve times. More... | |
| bool | is_doc_time_enabled () const |
| Is documentation of solve times enabled? More... | |
| bool | is_resolve_enabled () const |
| Boolean flag indicating if resolves are enabled. More... | |
| virtual void | enable_resolve () |
| virtual void | disable_resolve () |
| virtual void | solve (Problem *const &problem_pt, DoubleVector &result)=0 |
| virtual void | solve (DoubleMatrixBase *const &matrix_pt, const DoubleVector &rhs, DoubleVector &result) |
| virtual void | solve (DoubleMatrixBase *const &matrix_pt, const Vector< double > &rhs, Vector< double > &result) |
| virtual void | solve_transpose (Problem *const &problem_pt, DoubleVector &result) |
| virtual void | solve_transpose (DoubleMatrixBase *const &matrix_pt, const DoubleVector &rhs, DoubleVector &result) |
| virtual void | solve_transpose (DoubleMatrixBase *const &matrix_pt, const Vector< double > &rhs, Vector< double > &result) |
| virtual void | resolve (const DoubleVector &rhs, DoubleVector &result) |
| virtual void | resolve_transpose (const DoubleVector &rhs, DoubleVector &result) |
| virtual void | clean_up_memory () |
| virtual void | enable_computation_of_gradient () |
| void | disable_computation_of_gradient () |
| void | reset_gradient () |
| void | get_gradient (DoubleVector &gradient) |
| function to access the gradient, provided it has been computed More... | |
Public Member Functions inherited from oomph::DistributableLinearAlgebraObject | |
| DistributableLinearAlgebraObject () | |
| Default constructor - create a distribution. More... | |
| DistributableLinearAlgebraObject (const DistributableLinearAlgebraObject &matrix)=delete | |
| Broken copy constructor. More... | |
| void | operator= (const DistributableLinearAlgebraObject &)=delete |
| Broken assignment operator. More... | |
| virtual | ~DistributableLinearAlgebraObject () |
| Destructor. More... | |
| LinearAlgebraDistribution * | distribution_pt () const |
| access to the LinearAlgebraDistribution More... | |
| unsigned | nrow () const |
| access function to the number of global rows. More... | |
| unsigned | nrow_local () const |
| access function for the num of local rows on this processor. More... | |
| unsigned | nrow_local (const unsigned &p) const |
| access function for the num of local rows on this processor. More... | |
| unsigned | first_row () const |
| access function for the first row on this processor More... | |
| unsigned | first_row (const unsigned &p) const |
| access function for the first row on this processor More... | |
| bool | distributed () const |
| distribution is serial or distributed More... | |
| bool | distribution_built () const |
| void | build_distribution (const LinearAlgebraDistribution *const dist_pt) |
| void | build_distribution (const LinearAlgebraDistribution &dist) |
Static Protected Attributes | |
| static IdentityPreconditioner | Default_preconditioner |
Additional Inherited Members | |
Protected Member Functions inherited from oomph::DistributableLinearAlgebraObject | |
| void | clear_distribution () |
Base class for all linear iterative solvers. This merely defines standard interfaces for linear iterative solvers, so that different solvers can be used in a clean and transparent manner.
|
inline |
Constructor: Set (default) trivial preconditioner and set defaults for tolerance and max. number of iterations
References Default_preconditioner, Doc_convergence_history, First_time_solve_when_used_as_preconditioner, Max_iter, Preconditioner_pt, Setup_preconditioner_before_solve, Throw_error_after_max_iter, Tolerance, and Use_iterative_solver_as_preconditioner.
|
delete |
Broken copy constructor.
|
inlinevirtual |
|
inline |
Close convergence history output stream.
References Output_file_stream.
Referenced by PreconditionedFSICollapsibleChannelProblem< ELEMENT >::actions_after_newton_solve(), and PreconditionedFSICollapsibleChannelProblem< ELEMENT >::actions_before_newton_convergence_check().
|
inline |
Disable documentation of the convergence history.
References Doc_convergence_history.
|
inline |
Don't throw an error if we don't converge within max_iter (default).
References Throw_error_after_max_iter.
|
inline |
Disables the iterative solver be used as preconditioner (when calling the solve method it bypass the setup solver method – currently only used by Trilinos solver —)
References Use_iterative_solver_as_preconditioner.
|
inline |
Don't set up the preconditioner before the solve.
References Setup_preconditioner_before_solve.
|
inline |
Enable documentation of the convergence history.
References Doc_convergence_history.
Referenced by main().
|
inline |
Throw an error if we don't converge within max_iter.
References Throw_error_after_max_iter.
Referenced by oomph::MyProblem::build().
|
inline |
Enables the iterative solver be used as preconditioner (when calling the solve method it bypass the setup solver method – currently only used by Trilinos solver —)
References Use_iterative_solver_as_preconditioner.
|
inline |
Setup the preconditioner before the solve.
References Setup_preconditioner_before_solve.
|
pure virtual |
Number of iterations taken.
Implemented in oomph::GMRESBlockPreconditioner, oomph::GMRESBlockPreconditioner, oomph::HelmholtzGMRESMG< MATRIX >, oomph::ComplexGMRES< MATRIX >, oomph::ComplexDampedJacobi< MATRIX >, oomph::TrilinosAztecOOSolver, oomph::AugmentedProblemGMRES, oomph::GMRES< MATRIX >, oomph::DampedJacobi< MATRIX >, oomph::GS< CRDoubleMatrix >, oomph::GS< MATRIX >, oomph::BiCGStab< MATRIX >, oomph::CG< MATRIX >, and oomph::MGSolver< DIM >.
Referenced by oomph::MyProblem::actions_after_explicit_stage(), TiltedCavityProblem< ELEMENT >::actions_after_newton_step(), and oomph::MyProblem::actions_after_newton_step().
|
inlinevirtual |
returns the time taken to assemble the jacobian matrix and residual vector
Reimplemented from oomph::LinearSolver.
References Jacobian_setup_time.
|
inlinevirtual |
return the time taken to solve the linear system
Reimplemented from oomph::LinearSolver.
References Solution_time.
|
inline |
Access to max. number of iterations.
References Max_iter.
Referenced by oomph::MyProblem::build(), FSIChannelWithLeafletProblem< ELEMENT >::FSIChannelWithLeafletProblem(), main(), PreconditionedFSICollapsibleChannelProblem< ELEMENT >::PreconditionedFSICollapsibleChannelProblem(), run_it(), PMLStructuredCubicHelmholtz< ELEMENT >::set_gmres_multigrid_solver(), PMLHelmholtzMGProblem< ELEMENT >::set_gmres_multigrid_solver(), RefineableConvectionProblem< NST_ELEMENT, AD_ELEMENT >::switch_to_iterative_linear_solver(), TorusProblem< ELEMENT >::TorusProblem(), and TurekProblem< FLUID_ELEMENT, SOLID_ELEMENT >::TurekProblem().
|
inline |
Write convergence history into file with specified filename (automatically switches on doc). Optional second argument is a string that can be used (as a zone title) to identify what case we're running (e.g. what combination of linear solver and preconditioner or parameter values are used).
References Doc_convergence_history, Particles2023AnalysisHung::file_name, and Output_file_stream.
Referenced by PreconditionedFSICollapsibleChannelProblem< ELEMENT >::actions_before_newton_convergence_check(), MultiPoissonProblem< ELEMENT >::MultiPoissonProblem(), PMLStructuredCubicHelmholtz< ELEMENT >::set_gmres_multigrid_solver(), PMLHelmholtzMGProblem< ELEMENT >::set_gmres_multigrid_solver(), FluxPoissonMGProblem< ELEMENT, MESH >::set_multigrid_solver(), and UnitCubePoissonMGProblem< ELEMENT, MESH >::set_multigrid_solver().
|
delete |
Broken assignment operator.
|
inline |
Access function to preconditioner.
References Preconditioner_pt.
Referenced by oomph::AugmentedProblemGMRES::apply_schur_complement_preconditioner(), oomph::MyProblem::build(), oomph::TrilinosAztecOOSolver::clean_up_memory(), FSIChannelWithLeafletProblem< ELEMENT >::FSIChannelWithLeafletProblem(), main(), MultiPoissonProblem< ELEMENT >::MultiPoissonProblem(), parallel_test(), PreconditionedFSICollapsibleChannelProblem< ELEMENT >::PreconditionedFSICollapsibleChannelProblem(), PrescribedBoundaryDisplacementProblem< ELEMENT >::PrescribedBoundaryDisplacementProblem(), oomph::ProjectionProblem< PROJECTABLE_ELEMENT >::project(), PMLStructuredCubicHelmholtz< ELEMENT >::set_gmres_multigrid_solver(), PMLHelmholtzMGProblem< ELEMENT >::set_gmres_multigrid_solver(), FSIChannelWithLeafletProblem< ELEMENT >::set_iterative_solver(), PseudoElasticCollapsibleChannelProblem< FLUID_ELEMENT, SOLID_ELEMENT >::set_pseudo_elastic_fsi_solver(), set_pseudo_elastic_fsi_solver(), oomph::AugmentedProblemGMRES::solve_helper(), RefineableConvectionProblem< NST_ELEMENT, AD_ELEMENT >::switch_to_iterative_linear_solver(), TorusProblem< ELEMENT >::TorusProblem(), TurekProblem< FLUID_ELEMENT, SOLID_ELEMENT >::TurekProblem(), oomph::GMRES< MATRIX >::update(), oomph::HelmholtzGMRESMG< MATRIX >::update(), and RefineableConvectionProblem< NST_ELEMENT, AD_ELEMENT >::~RefineableConvectionProblem().
|
inline |
Access function to preconditioner (const version)
References Preconditioner_pt.
|
inlinevirtual |
returns the the time taken to setup the preconditioner
Reimplemented in oomph::GS< CRDoubleMatrix >, and oomph::GS< MATRIX >.
References Preconditioner_setup_time.
Referenced by oomph::MyProblem::actions_after_explicit_stage(), and oomph::MyProblem::actions_after_newton_step().
|
inline |
Access to convergence tolerance.
References Tolerance.
Referenced by oomph::ComplexDampedJacobi< MATRIX >::calculate_omega(), oomph::ComplexGMRES< MATRIX >::generate_plane_rotation(), oomph::HelmholtzGMRESMG< MATRIX >::generate_plane_rotation(), main(), MultiPoissonProblem< ELEMENT >::MultiPoissonProblem(), PreconditionedFSICollapsibleChannelProblem< ELEMENT >::PreconditionedFSICollapsibleChannelProblem(), PMLStructuredCubicHelmholtz< ELEMENT >::set_gmres_multigrid_solver(), PMLHelmholtzMGProblem< ELEMENT >::set_gmres_multigrid_solver(), set_pseudo_elastic_fsi_solver(), oomph::MGSolver< DIM >::setup_smoothers(), oomph::HelmholtzMGPreconditioner< DIM >::setup_smoothers(), RefineableConvectionProblem< NST_ELEMENT, AD_ELEMENT >::switch_to_iterative_linear_solver(), TiltedCavityProblem< ELEMENT >::TiltedCavityProblem(), TorusProblem< ELEMENT >::TorusProblem(), and TurekProblem< FLUID_ELEMENT, SOLID_ELEMENT >::TurekProblem().
|
staticprotected |
Default preconditioner: The base class for preconditioners is a fully functional (if trivial!) preconditioner.
Default preconditioner for iterative solvers: The base class for preconditioners is a fully functional (if trivial!) preconditioner.
Referenced by IterativeLinearSolver().
|
protected |
Flag indicating if the convergence history is to be documented
Referenced by disable_doc_convergence_history(), enable_doc_convergence_history(), IterativeLinearSolver(), open_convergence_history_file_stream(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), and oomph::AugmentedProblemGMRES::solve_helper().
|
protected |
When the iterative solver is used a preconditioner then we call the setup of solver method only once (the first time the solve method is called)
Referenced by IterativeLinearSolver(), and oomph::TrilinosAztecOOSolver::solve().
|
protected |
Jacobian setup time.
Referenced by jacobian_setup_time(), oomph::GS< CRDoubleMatrix >::solve(), oomph::AugmentedProblemGMRES::solve(), oomph::HelmholtzGMRESMG< MATRIX >::solve(), and oomph::HelmholtzFGMRESMG< MATRIX >::solve().
|
protected |
Maximum number of iterations.
Referenced by IterativeLinearSolver(), max_iter(), oomph::TrilinosAztecOOSolver::max_iter(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), oomph::AugmentedProblemGMRES::solve_helper(), oomph::TrilinosAztecOOSolver::solve_using_AztecOO(), and oomph::TrilinosAztecOOSolver::solver_setup().
|
protected |
Output file stream for convergence history.
Referenced by close_convergence_history_file_stream(), open_convergence_history_file_stream(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), and oomph::AugmentedProblemGMRES::solve_helper().
|
protected |
Pointer to the preconditioner.
Referenced by IterativeLinearSolver(), preconditioner_pt(), and oomph::TrilinosAztecOOSolver::solver_setup().
|
protected |
Preconditioner setup time.
Referenced by preconditioner_setup_time(), and oomph::AugmentedProblemGMRES::solve_helper().
|
protected |
indicates whether the preconditioner should be setup before solve. Default = true;
Referenced by disable_setup_preconditioner_before_solve(), enable_setup_preconditioner_before_solve(), IterativeLinearSolver(), oomph::AugmentedProblemGMRES::solve_helper(), and oomph::TrilinosAztecOOSolver::solver_setup().
|
protected |
linear solver solution time
Referenced by linear_solver_solution_time(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), and oomph::AugmentedProblemGMRES::solve_helper().
|
protected |
Should we throw an error instead of just returning when we hit the max iterations?
Referenced by disable_error_after_max_iter(), enable_error_after_max_iter(), IterativeLinearSolver(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), and oomph::AugmentedProblemGMRES::solve_helper().
|
protected |
Convergence tolerance.
Referenced by IterativeLinearSolver(), oomph::MGSolver< DIM >::MGSolver(), oomph::GMRESBlockPreconditioner::preconditioner_solve(), oomph::GS< CRDoubleMatrix >::solve_helper(), oomph::AugmentedProblemGMRES::solve_helper(), oomph::TrilinosAztecOOSolver::solve_using_AztecOO(), tolerance(), oomph::TrilinosAztecOOSolver::tolerance(), and oomph::TrilinosAztecOOSolver::TrilinosAztecOOSolver().
|
protected |
Use the iterative solver as preconditioner.
Referenced by disable_iterative_solver_as_preconditioner(), enable_iterative_solver_as_preconditioner(), IterativeLinearSolver(), and oomph::TrilinosAztecOOSolver::solve().