Anasazi Version of the Day
Loading...
Searching...
No Matches
AnasaziTraceMinSolMgr.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Anasazi: Block Eigensolvers Package
4//
5// Copyright 2004 NTESS and the Anasazi contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10#ifndef ANASAZI_TRACEMIN_SOLMGR_HPP
11#define ANASAZI_TRACEMIN_SOLMGR_HPP
12
19#include "AnasaziConfigDefs.hpp"
20#include "AnasaziTypes.hpp"
21
24
25#include "AnasaziTraceMin.hpp"
27#include "AnasaziBasicSort.hpp"
33#include "Teuchos_BLAS.hpp"
34#include "Teuchos_LAPACK.hpp"
35#include "Teuchos_TimeMonitor.hpp"
36#ifdef TEUCHOS_DEBUG
37# include <Teuchos_FancyOStream.hpp>
38#endif
39
40namespace Anasazi {
41namespace Experimental {
42
43template<class ScalarType, class MV, class OP>
44
77class TraceMinSolMgr : public TraceMinBaseSolMgr<ScalarType,MV,OP> {
78
79 private:
82 typedef Teuchos::ScalarTraits<ScalarType> SCT;
83 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
84 typedef Teuchos::ScalarTraits<MagnitudeType> MT;
85
86 public:
87
89
90
101 TraceMinSolMgr( const Teuchos::RCP<Eigenproblem<ScalarType,MV,OP> > &problem,
102 Teuchos::ParameterList &pl );
104
105 private:
106
107 int maxits_;
108
109 // Test whether we have exceeded the maximum number of iterations
110 bool exceededMaxIter() { return (this->iter_ >= maxits_); };
111
112 // TraceMin does not restart, so this will always return false
113 bool needToRestart(const Teuchos::RCP< TraceMinBase<ScalarType,MV,OP> > solver) { return false; };
114
115 // TraceMin does not restart, so this will throw an exception
116 bool performRestart(int &numRestarts, Teuchos::RCP< TraceMinBase<ScalarType,MV,OP> > solver)
117 { TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Anasazi::TraceMinSolMgr::performRestart(): TraceMin does not perform restarts!"); };
118
119 // Returns a new TraceMin solver object
120 Teuchos::RCP< TraceMinBase<ScalarType,MV,OP> > createSolver(
121 const Teuchos::RCP<SortManager<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> > &sorter,
122 const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &outputtest,
123 const Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> > &ortho,
124 Teuchos::ParameterList &plist
125 );
126};
127
128
129//---------------------------------------------------------------------------//
130// Prevent instantiation on complex scalar type
131// FIXME: this really is just a current flaw in the implementation, TraceMin
132// *should* work for Hermitian matrices
133//---------------------------------------------------------------------------//
134template <class MagnitudeType, class MV, class OP>
135class TraceMinSolMgr<std::complex<MagnitudeType>,MV,OP>
136{
137 public:
138
139 typedef std::complex<MagnitudeType> ScalarType;
141 const RCP<Eigenproblem<ScalarType,MV,OP> > &problem,
142 Teuchos::ParameterList &pl )
143 {
144 // Provide a compile error when attempting to instantiate on complex type
145 MagnitudeType::this_class_is_missing_a_specialization();
146 }
147};
148
150// Constructor - accepts maximum iterations in addition to the other parameters of the abstract base class
151template<class ScalarType, class MV, class OP>
152TraceMinSolMgr<ScalarType,MV,OP>::TraceMinSolMgr( const Teuchos::RCP<Eigenproblem<ScalarType,MV,OP> > &problem, Teuchos::ParameterList &pl ) :
153 TraceMinBaseSolMgr<ScalarType,MV,OP>(problem,pl)
154{
155 // Get the maximum number of iterations
156 maxits_ = pl.get("Maximum Iterations", 100);
157 TEUCHOS_TEST_FOR_EXCEPTION(maxits_ < 1, std::invalid_argument, "Anasazi::TraceMinSolMgr::constructor(): \"Maximum Iterations\" must be strictly positive.");
158
159 // block size: default is 2* nev()
160 // TODO: Find out minimum value
161 this->blockSize_ = pl.get("Block Size",2*this->problem_->getNEV());
162 TEUCHOS_TEST_FOR_EXCEPTION(this->blockSize_ < this->problem_->getNEV(), std::invalid_argument,
163 "Anasazi::TraceMinSolMgr::constructor(): \"Block Size\" must be greater than or equal to the number of desired eigenpairs.");
164
165 this->useHarmonic_ = pl.get("Use Harmonic Ritz Values", false);
166 TEUCHOS_TEST_FOR_EXCEPTION(this->useHarmonic_, std::invalid_argument,
167 "Anasazi::TraceMinSolMgr::constructor(): Please disable the harmonic Ritz values. It doesn't make sense to use them with TraceMin, which does not use expanding subspaces. Perhaps you wanted TraceMin-Davidson?");
168
169 // TraceMin does not restart, so the number of blocks and number of restart blocks will always be 1
170 this->numBlocks_ = 1;
171 this->numRestartBlocks_ = 1;
172
173 TEUCHOS_TEST_FOR_EXCEPTION(static_cast<ptrdiff_t>(this->numBlocks_)*this->blockSize_ + this->maxLocked_ > MVT::GetGlobalLength(*this->problem_->getInitVec()),
174 std::invalid_argument,
175 "Anasazi::TraceMinSolMgr::constructor(): Potentially impossible orthogonality requests. Reduce basis size or locking size.");
176
177 TEUCHOS_TEST_FOR_EXCEPTION(this->maxLocked_ + this->blockSize_ < this->problem_->getNEV(), std::invalid_argument,
178 "Anasazi::TraceMinDavidsonSolMgr: Not enough storage space for requested number of eigenpairs.");
179}
180
181
183// Returns a new TraceMin solver object
184template <class ScalarType, class MV, class OP>
185Teuchos::RCP< TraceMinBase<ScalarType,MV,OP> > TraceMinSolMgr<ScalarType,MV,OP>::createSolver(
186 const Teuchos::RCP<SortManager<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> > &sorter,
187 const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &outputtest,
188 const Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> > &ortho,
189 Teuchos::ParameterList &plist
190 )
191{
192 return Teuchos::rcp( new TraceMin<ScalarType,MV,OP>(this->problem_,sorter,this->printer_,outputtest,ortho,plist) );
193}
194
195
196}} // end Anasazi namespace
197
198#endif /* ANASAZI_TRACEMIN_SOLMGR_HPP */
Basic implementation of the Anasazi::SortManager class.
Anasazi header file which uses auto-configuration information to include necessary C++ headers.
Abstract base class which defines the interface required by an eigensolver and status test class to c...
Abstract class definition for Anasazi Output Managers.
Class which provides internal utilities for the Anasazi solvers.
Status test for forming logical combinations of other status tests.
Special StatusTest for printing status tests.
A status test for testing the norm of the eigenvectors residuals.
A status test for testing the norm of the eigenvectors residuals along with a set of auxiliary eigenv...
The Anasazi::TraceMinBaseSolMgr provides an abstract base class for the TraceMin series of solver man...
Implementation of the trace minimization eigensolver.
Types and exceptions used within Anasazi solvers and interfaces.
This class defines the interface required by an eigensolver and status test class to compute solution...
The Anasazi::TraceMinBaseSolMgr provides an abstract base class for the TraceMin series of solver man...
This is an abstract base class for the trace minimization eigensolvers.
The Anasazi::TraceMinSolMgr provides a flexible solver manager over the TraceMin eigensolver.
TraceMinSolMgr(const Teuchos::RCP< Eigenproblem< ScalarType, MV, OP > > &problem, Teuchos::ParameterList &pl)
Basic constructor for TraceMinSolMgr.
This class implements a TraceMIN iteration, a preconditioned iteration for solving linear symmetric p...
Anasazi's templated virtual class for providing routines for orthogonalization and orthonormalization...
Traits class which defines basic operations on multivectors.
Virtual base class which defines basic traits for the operator type.
Anasazi's templated pure virtual class for managing the sorting of approximate eigenvalues computed b...
Common interface of stopping criteria for Anasazi's solvers.
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package.