PyTrilinos::Anasazi::LOBPCGEpetra Class Reference

Inheritance diagram for PyTrilinos::Anasazi::LOBPCGEpetra:

Inheritance graph
[legend]
Collaboration diagram for PyTrilinos::Anasazi::LOBPCGEpetra:

Collaboration graph
[legend]

List of all members.

Public Member Functions

def __init__
def iterate
def initialize
def isInitialized
def getState
def getNumIters
def resetNumIters
def getRitzVectors
def getRitzValues
def getRitzIndex
def getResNorms
def getRes2Norms
def getRitzRes2Norms
def getCurSubspaceDim
def getMaxSubspaceDim
def setStatusTest
def getStatusTest
def getProblem
def setBlockSize
def getBlockSize
def setAuxVecs
def getAuxVecs
def setFullOrtho
def getFullOrtho
def hasP
def currentStatus
def __init__
def iterate
def initialize
def isInitialized
def getState
def getNumIters
def resetNumIters
def getRitzVectors
def getRitzValues
def getRitzIndex
def getResNorms
def getRes2Norms
def getRitzRes2Norms
def getCurSubspaceDim
def getMaxSubspaceDim
def setStatusTest
def getStatusTest
def getProblem
def setBlockSize
def getBlockSize
def setAuxVecs
def getAuxVecs
def setFullOrtho
def getFullOrtho
def hasP
def currentStatus

Public Attributes

 this


Detailed Description

This class provides the Locally Optimal Block Preconditioned Conjugate
Gradient (LOBPCG) iteration, a preconditioned iteration for solving
linear Hermitian eigenproblems.

This implementation is a modification of the one found in A. Knyazev,
"Toward the optimal preconditioned eigensolver: Locally optimal block
preconditioner conjugate gradient method", SIAM J. Sci. Comput., vol
23, n 2, pp. 517-541.

The modification consists of the orthogonalization steps recommended
in U. Hetmaniuk and R. Lehoucq, "Basis Selection in LOBPCG", Journal
of Computational Physics.

These modifcation are referred to as full orthogonalization, and
consist of also conducting the local optimization using an orthonormal
basis.

Chris Baker, Ulrich Hetmaniuk, Rich Lehoucq, Heidi Thornquist

C++ includes: AnasaziLOBPCG.hpp 

Member Function Documentation

def PyTrilinos::Anasazi::LOBPCGEpetra::__init__ (   self,
  args 
)

__init__(self, Teuchos::RCP<(Anasazi::Eigenproblem<(double,Epetra_MultiVector,Epetra_Operator)>)> problem, 
    Teuchos::RCP<(Anasazi::SortManager<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)> sorter, 
    Teuchos::RCP<(Anasazi::OutputManager<(double)>)> printer, 
    Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)> tester, 
    Teuchos::RCP<(Anasazi::MatOrthoManager<(double,Epetra_MultiVector,Epetra_Operator)>)> ortho, 
    ParameterList params) -> LOBPCGEpetra

Anasazi::LOBPCG<
ScalarType, MV, OP >::LOBPCG(const Teuchos::RCP< Eigenproblem<
ScalarType, MV, OP > > &problem, const Teuchos::RCP< SortManager<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > >
&sorter, const Teuchos::RCP< OutputManager< ScalarType > > &printer,
const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, const
Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > &ortho,
Teuchos::ParameterList &params)

LOBPCG constructor with eigenproblem, solver utilities, and parameter
list of solver options.

This constructor takes pointers required by the eigensolver, in
addition to a parameter list of options for the eigensolver. These
options include the following: "Block Size" - an int specifying the
block size used by the algorithm. This can also be specified using the
setBlockSize() method.

"Full Ortho" - a bool specifying whether the solver should employ a
full orthogonalization technique. This can also be specified using the
setFullOrtho() method. 

def PyTrilinos::Anasazi::LOBPCGEpetra::__init__ (   self,
  args 
)

__init__(self, Teuchos::RCP<(Anasazi::Eigenproblem<(double,Epetra_MultiVector,Epetra_Operator)>)> problem, 
    Teuchos::RCP<(Anasazi::SortManager<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)> sorter, 
    Teuchos::RCP<(Anasazi::OutputManager<(double)>)> printer, 
    Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)> tester, 
    Teuchos::RCP<(Anasazi::MatOrthoManager<(double,Epetra_MultiVector,Epetra_Operator)>)> ortho, 
    ParameterList params) -> LOBPCGEpetra

Anasazi::LOBPCG<
ScalarType, MV, OP >::LOBPCG(const Teuchos::RCP< Eigenproblem<
ScalarType, MV, OP > > &problem, const Teuchos::RCP< SortManager<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > >
&sorter, const Teuchos::RCP< OutputManager< ScalarType > > &printer,
const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, const
Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > &ortho,
Teuchos::ParameterList &params)

LOBPCG constructor with eigenproblem, solver utilities, and parameter
list of solver options.

This constructor takes pointers required by the eigensolver, in
addition to a parameter list of options for the eigensolver. These
options include the following: "Block Size" - an int specifying the
block size used by the algorithm. This can also be specified using the
setBlockSize() method.

"Full Ortho" - a bool specifying whether the solver should employ a
full orthogonalization technique. This can also be specified using the
setFullOrtho() method. 

def PyTrilinos::Anasazi::LOBPCGEpetra::currentStatus (   self,
  args 
)

currentStatus(self, ostream os)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::currentStatus(std::ostream &os)

This method requests that the solver print out its current status to
screen. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::currentStatus (   self,
  args 
)

currentStatus(self, ostream os)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::currentStatus(std::ostream &os)

This method requests that the solver print out its current status to
screen. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getAuxVecs (   self,
  args 
)

getAuxVecs(self) -> Teuchos::Array<(Teuchos::RCP<(q(const).Epetra_MultiVector)>)>

Teuchos::Array<
Teuchos::RCP< const MV > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getAuxVecs() const

Get the current auxiliary vectors. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getAuxVecs (   self,
  args 
)

getAuxVecs(self) -> Teuchos::Array<(Teuchos::RCP<(q(const).Epetra_MultiVector)>)>

Teuchos::Array<
Teuchos::RCP< const MV > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getAuxVecs() const

Get the current auxiliary vectors. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getBlockSize (   self,
  args 
)

getBlockSize(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getBlockSize() const

Get the blocksize to be used by the iterative solver in solving this
eigenproblem. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getBlockSize (   self,
  args 
)

getBlockSize(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getBlockSize() const

Get the blocksize to be used by the iterative solver in solving this
eigenproblem. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getCurSubspaceDim (   self,
  args 
)

getCurSubspaceDim(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getCurSubspaceDim() const

Get the dimension of the search subspace used to generate the current
eigenvectors and eigenvalues.

LOBPCG employs a sequential subspace iteration, maintaining a fixed-
rank basis, as opposed to an expanding subspace mechanism employed by
Krylov-subspace solvers like BlockKrylovSchur and BlockDavidson.

An integer specifying the rank of the subspace generated by the
eigensolver. If isInitialized() == false, the return is 0. Otherwise,
the return will be 2*getBlockSize() or 3*getBlockSize(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getCurSubspaceDim (   self,
  args 
)

getCurSubspaceDim(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getCurSubspaceDim() const

Get the dimension of the search subspace used to generate the current
eigenvectors and eigenvalues.

LOBPCG employs a sequential subspace iteration, maintaining a fixed-
rank basis, as opposed to an expanding subspace mechanism employed by
Krylov-subspace solvers like BlockKrylovSchur and BlockDavidson.

An integer specifying the rank of the subspace generated by the
eigensolver. If isInitialized() == false, the return is 0. Otherwise,
the return will be 2*getBlockSize() or 3*getBlockSize(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getFullOrtho (   self,
  args 
)

getFullOrtho(self) -> bool

bool
Anasazi::LOBPCG< ScalarType, MV, OP >::getFullOrtho() const

Determine if the LOBPCG iteration is using full orthogonality. 

def PyTrilinos::Anasazi::LOBPCGEpetra::getFullOrtho (   self,
  args 
)

getFullOrtho(self) -> bool

bool
Anasazi::LOBPCG< ScalarType, MV, OP >::getFullOrtho() const

Determine if the LOBPCG iteration is using full orthogonality. 

def PyTrilinos::Anasazi::LOBPCGEpetra::getMaxSubspaceDim (   self,
  args 
)

getMaxSubspaceDim(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getMaxSubspaceDim() const

Get the maximum dimension allocated for the search subspace. For
LOBPCG, this always returns 3*getBlockSize(), the dimension of the
subspace colspan([X H P]). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getMaxSubspaceDim (   self,
  args 
)

getMaxSubspaceDim(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getMaxSubspaceDim() const

Get the maximum dimension allocated for the search subspace. For
LOBPCG, this always returns 3*getBlockSize(), the dimension of the
subspace colspan([X H P]). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getNumIters (   self,
  args 
)

getNumIters(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getNumIters() const

Get the current iteration count. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getNumIters (   self,
  args 
)

getNumIters(self) -> int

int
Anasazi::LOBPCG< ScalarType, MV, OP >::getNumIters() const

Get the current iteration count. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getProblem (   self,
  args 
)

getProblem(self) -> EigenproblemEpetra

const
Eigenproblem< ScalarType, MV, OP > & Anasazi::LOBPCG< ScalarType, MV,
OP >::getProblem() const

Get a constant reference to the eigenvalue problem. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getProblem (   self,
  args 
)

getProblem(self) -> EigenproblemEpetra

const
Eigenproblem< ScalarType, MV, OP > & Anasazi::LOBPCG< ScalarType, MV,
OP >::getProblem() const

Get a constant reference to the eigenvalue problem. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRes2Norms (   self,
  args 
)

getRes2Norms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType >
Anasazi::LOBPCG< ScalarType, MV, OP >::getRes2Norms()

Get the current residual 2-norms.

A vector of length getBlockSize() containing the 2-norms of the
residuals. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRes2Norms (   self,
  args 
)

getRes2Norms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType >
Anasazi::LOBPCG< ScalarType, MV, OP >::getRes2Norms()

Get the current residual 2-norms.

A vector of length getBlockSize() containing the 2-norms of the
residuals. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getResNorms (   self,
  args 
)

getResNorms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType >
Anasazi::LOBPCG< ScalarType, MV, OP >::getResNorms()

Get the current residual norms.

A vector of length getBlockSize() containing the norms of the
residuals, with respect to the orthogonalization manager norm()
method. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getResNorms (   self,
  args 
)

getResNorms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector<
typename Teuchos::ScalarTraits< ScalarType >::magnitudeType >
Anasazi::LOBPCG< ScalarType, MV, OP >::getResNorms()

Get the current residual norms.

A vector of length getBlockSize() containing the norms of the
residuals, with respect to the orthogonalization manager norm()
method. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzIndex (   self,
  args 
)

getRitzIndex(self) -> VectorInt

std::vector<
int > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzIndex()

Get the index used for extracting Ritz vectors from getRitzVectors().

Because BlockDavidson is a Hermitian solver, all Ritz values are real
and all Ritz vectors can be represented in a single column of a
multivector. Therefore, getRitzIndex() is not needed when using the
output from getRitzVectors().

An int vector of size getCurSubspaceDim() composed of zeros. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzIndex (   self,
  args 
)

getRitzIndex(self) -> VectorInt

std::vector<
int > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzIndex()

Get the index used for extracting Ritz vectors from getRitzVectors().

Because BlockDavidson is a Hermitian solver, all Ritz values are real
and all Ritz vectors can be represented in a single column of a
multivector. Therefore, getRitzIndex() is not needed when using the
output from getRitzVectors().

An int vector of size getCurSubspaceDim() composed of zeros. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzRes2Norms (   self,
  args 
)

getRitzRes2Norms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector< typename Teuchos::ScalarTraits< ScalarType
>::magnitudeType > Anasazi::LOBPCG< ScalarType, MV, OP
>::getRitzRes2Norms()

Get the 2-norms of the residuals.

The Ritz residuals are not defined for the LOBPCG iteration. Hence,
this method returns the 2-norms of the direct residuals, and is
equivalent to calling getRes2Norms().

A vector of length getBlockSize() containing the 2-norms of the direct
residuals. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzRes2Norms (   self,
  args 
)

getRitzRes2Norms(self) -> std::vector<(Teuchos::ScalarTraits<(double)>::magnitudeType,std::allocator<(Teuchos::ScalarTraits<(double)>::magnitudeType)>)>

std::vector< typename Teuchos::ScalarTraits< ScalarType
>::magnitudeType > Anasazi::LOBPCG< ScalarType, MV, OP
>::getRitzRes2Norms()

Get the 2-norms of the residuals.

The Ritz residuals are not defined for the LOBPCG iteration. Hence,
this method returns the 2-norms of the direct residuals, and is
equivalent to calling getRes2Norms().

A vector of length getBlockSize() containing the 2-norms of the direct
residuals. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzValues (   self,
  args 
)

getRitzValues(self) -> VectorValueDouble

std::vector<
Value< ScalarType > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getRitzValues()

Get the Ritz values from the previous iteration.

A vector of length getCurSubspaceDim() containing the Ritz values from
the previous projected eigensolve. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzValues (   self,
  args 
)

getRitzValues(self) -> VectorValueDouble

std::vector<
Value< ScalarType > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getRitzValues()

Get the Ritz values from the previous iteration.

A vector of length getCurSubspaceDim() containing the Ritz values from
the previous projected eigensolve. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzVectors (   self,
  args 
)

getRitzVectors(self) -> Teuchos::RCP<(q(const).Epetra_MultiVector)>

Teuchos::RCP<
const MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzVectors()

Get the Ritz vectors from the previous iteration.

A multivector with getBlockSize() vectors containing the sorted Ritz
vectors corresponding to the most significant Ritz values. The i-th
vector of the return corresponds to the i-th Ritz vector; there is no
need to use getRitzIndex(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getRitzVectors (   self,
  args 
)

getRitzVectors(self) -> Teuchos::RCP<(q(const).Epetra_MultiVector)>

Teuchos::RCP<
const MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzVectors()

Get the Ritz vectors from the previous iteration.

A multivector with getBlockSize() vectors containing the sorted Ritz
vectors corresponding to the most significant Ritz values. The i-th
vector of the return corresponds to the i-th Ritz vector; there is no
need to use getRitzIndex(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getState (   self,
  args 
)

getState(self) -> Anasazi::LOBPCGState<(double,Epetra_MultiVector)>

LOBPCGState<
ScalarType, MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getState()
const

Get the current state of the eigensolver.

The data is only valid if isInitialized() == true. The data for the
search directions P is only meaningful if hasP() == true. Finally, the
data for the preconditioned residual (H) is only meaningful in the
situation where the solver throws an LOBPCGRitzFailure exception
during iterate().

An LOBPCGState object containing const views to the current solver
state. 

def PyTrilinos::Anasazi::LOBPCGEpetra::getState (   self,
  args 
)

getState(self) -> Anasazi::LOBPCGState<(double,Epetra_MultiVector)>

LOBPCGState<
ScalarType, MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getState()
const

Get the current state of the eigensolver.

The data is only valid if isInitialized() == true. The data for the
search directions P is only meaningful if hasP() == true. Finally, the
data for the preconditioned residual (H) is only meaningful in the
situation where the solver throws an LOBPCGRitzFailure exception
during iterate().

An LOBPCGState object containing const views to the current solver
state. 

def PyTrilinos::Anasazi::LOBPCGEpetra::getStatusTest (   self,
  args 
)

getStatusTest(self) -> Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)>

Teuchos::RCP<
StatusTest< ScalarType, MV, OP > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getStatusTest() const

Get the current StatusTest used by the solver. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::getStatusTest (   self,
  args 
)

getStatusTest(self) -> Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)>

Teuchos::RCP<
StatusTest< ScalarType, MV, OP > > Anasazi::LOBPCG< ScalarType, MV, OP
>::getStatusTest() const

Get the current StatusTest used by the solver. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::hasP (   self,
  args 
)

hasP(self) -> bool

bool Anasazi::LOBPCG<
ScalarType, MV, OP >::hasP()

Indicates whether the search direction given by getState() is valid.

def PyTrilinos::Anasazi::LOBPCGEpetra::hasP (   self,
  args 
)

hasP(self) -> bool

bool Anasazi::LOBPCG<
ScalarType, MV, OP >::hasP()

Indicates whether the search direction given by getState() is valid.

def PyTrilinos::Anasazi::LOBPCGEpetra::initialize (   self,
  args 
)

initialize(self, Anasazi::LOBPCGState<(double,Epetra_MultiVector)> newstate)
initialize(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::initialize()

Initialize the solver with the initial vectors from the eigenproblem
or random data. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::initialize (   self,
  args 
)

initialize(self, Anasazi::LOBPCGState<(double,Epetra_MultiVector)> newstate)
initialize(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::initialize()

Initialize the solver with the initial vectors from the eigenproblem
or random data. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::isInitialized (   self,
  args 
)

isInitialized(self) -> bool

bool
Anasazi::LOBPCG< ScalarType, MV, OP >::isInitialized() const

Indicates whether the solver has been initialized or not.

bool indicating the state of the solver.

If isInitialized() == true: X is orthogonal to auxiliary vectors and
has orthonormal columns

KX == Op*X

MX == M*X if M != Teuchos::null  Otherwise, MX == Teuchos::null

getRitzValues() returns the sorted Ritz values with respect to X

getResNorms(), getRes2Norms(), getRitzResNorms() are correct

If hasP() == true, P orthogonal to auxiliary vectors

If getFullOrtho() == true, P is orthogonal to X and has orthonormal
columns

KP == Op*P

MP == M*P if M != Teuchos::null  Otherwise, MP == Teuchos::null 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::isInitialized (   self,
  args 
)

isInitialized(self) -> bool

bool
Anasazi::LOBPCG< ScalarType, MV, OP >::isInitialized() const

Indicates whether the solver has been initialized or not.

bool indicating the state of the solver.

If isInitialized() == true: X is orthogonal to auxiliary vectors and
has orthonormal columns

KX == Op*X

MX == M*X if M != Teuchos::null  Otherwise, MX == Teuchos::null

getRitzValues() returns the sorted Ritz values with respect to X

getResNorms(), getRes2Norms(), getRitzResNorms() are correct

If hasP() == true, P orthogonal to auxiliary vectors

If getFullOrtho() == true, P is orthogonal to X and has orthonormal
columns

KP == Op*P

MP == M*P if M != Teuchos::null  Otherwise, MP == Teuchos::null 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::iterate (   self,
  args 
)

iterate(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::iterate()

This method performs LOBPCG iterations until the status test indicates
the need to stop or an error occurs (in which case, an exception is
thrown).

iterate() will first determine whether the solver is initialized; if
not, it will call initialize() using default arguments. After
initialization, the solver performs LOBPCG iterations until the status
test evaluates as Passed, at which point the method returns to the
caller.

The LOBPCG iteration proceeds as follows: The current residual (R) is
preconditioned to form H

H is orthogonalized against the auxiliary vectors and, if full
orthogonalization  is enabled, against X and P.

The basis [X H P] is used to project the problem matrices.

The projected eigenproblem is solved, and the desired eigenvectors and
eigenvalues are selected.

These are used to form the new eigenvector estimates (X) and the
search directions (P).  If full orthogonalization is enabled, these
are generated to be mutually orthogonal and with orthonormal columns.

The new residual (R) is formed.

The status test is queried at the beginning of the iteration.

Possible exceptions thrown include std::logic_error,
std::invalid_argument or one of the LOBPCG-specific exceptions. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::iterate (   self,
  args 
)

iterate(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::iterate()

This method performs LOBPCG iterations until the status test indicates
the need to stop or an error occurs (in which case, an exception is
thrown).

iterate() will first determine whether the solver is initialized; if
not, it will call initialize() using default arguments. After
initialization, the solver performs LOBPCG iterations until the status
test evaluates as Passed, at which point the method returns to the
caller.

The LOBPCG iteration proceeds as follows: The current residual (R) is
preconditioned to form H

H is orthogonalized against the auxiliary vectors and, if full
orthogonalization  is enabled, against X and P.

The basis [X H P] is used to project the problem matrices.

The projected eigenproblem is solved, and the desired eigenvectors and
eigenvalues are selected.

These are used to form the new eigenvector estimates (X) and the
search directions (P).  If full orthogonalization is enabled, these
are generated to be mutually orthogonal and with orthonormal columns.

The new residual (R) is formed.

The status test is queried at the beginning of the iteration.

Possible exceptions thrown include std::logic_error,
std::invalid_argument or one of the LOBPCG-specific exceptions. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::resetNumIters (   self,
  args 
)

resetNumIters(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::resetNumIters()

Reset the iteration count. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::resetNumIters (   self,
  args 
)

resetNumIters(self)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::resetNumIters()

Reset the iteration count. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setAuxVecs (   self,
  args 
)

setAuxVecs(self, Teuchos::Array<(Teuchos::RCP<(q(const).Epetra_MultiVector)>)> auxvecs)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setAuxVecs(const
Teuchos::Array< Teuchos::RCP< const MV > > &auxvecs)

Set the auxiliary vectors for the solver.

Because the current iterate X and search direction P cannot be assumed
orthogonal to the new auxiliary vectors, a call to setAuxVecs() with a
non-empty argument will reset the solver to the uninitialized state.

In order to preserve the current state, the user will need to extract
it from the solver using getState(), orthogonalize it against the new
auxiliary vectors, and manually reinitialize the solver using
initialize(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setAuxVecs (   self,
  args 
)

setAuxVecs(self, Teuchos::Array<(Teuchos::RCP<(q(const).Epetra_MultiVector)>)> auxvecs)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setAuxVecs(const
Teuchos::Array< Teuchos::RCP< const MV > > &auxvecs)

Set the auxiliary vectors for the solver.

Because the current iterate X and search direction P cannot be assumed
orthogonal to the new auxiliary vectors, a call to setAuxVecs() with a
non-empty argument will reset the solver to the uninitialized state.

In order to preserve the current state, the user will need to extract
it from the solver using getState(), orthogonalize it against the new
auxiliary vectors, and manually reinitialize the solver using
initialize(). 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setBlockSize (   self,
  args 
)

setBlockSize(self, int blockSize)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setBlockSize(int blockSize)

Set the blocksize to be used by the iterative solver in solving this
eigenproblem.

If the block size is reduced, then the new iterate (and residual and
search direction) are chosen as the subset of the current iterate
preferred by the sort manager. Otherwise, the solver state is set to
uninitialized. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setBlockSize (   self,
  args 
)

setBlockSize(self, int blockSize)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setBlockSize(int blockSize)

Set the blocksize to be used by the iterative solver in solving this
eigenproblem.

If the block size is reduced, then the new iterate (and residual and
search direction) are chosen as the subset of the current iterate
preferred by the sort manager. Otherwise, the solver state is set to
uninitialized. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setFullOrtho (   self,
  args 
)

setFullOrtho(self, bool fullOrtho)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setFullOrtho(bool fullOrtho)

Instruct the LOBPCG iteration to use full orthogonality.

If the getFullOrtho() == false and isInitialized() == true and hasP()
== true, then P will be invalidated by setting full orthogonalization
to true. 

def PyTrilinos::Anasazi::LOBPCGEpetra::setFullOrtho (   self,
  args 
)

setFullOrtho(self, bool fullOrtho)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setFullOrtho(bool fullOrtho)

Instruct the LOBPCG iteration to use full orthogonality.

If the getFullOrtho() == false and isInitialized() == true and hasP()
== true, then P will be invalidated by setting full orthogonalization
to true. 

def PyTrilinos::Anasazi::LOBPCGEpetra::setStatusTest (   self,
  args 
)

setStatusTest(self, Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)> test)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setStatusTest(Teuchos::RCP<
StatusTest< ScalarType, MV, OP > > test)

Set a new StatusTest for the solver. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.

def PyTrilinos::Anasazi::LOBPCGEpetra::setStatusTest (   self,
  args 
)

setStatusTest(self, Teuchos::RCP<(Anasazi::StatusTest<(double,Epetra_MultiVector,Epetra_Operator)>)> test)

void
Anasazi::LOBPCG< ScalarType, MV, OP >::setStatusTest(Teuchos::RCP<
StatusTest< ScalarType, MV, OP > > test)

Set a new StatusTest for the solver. 

Reimplemented from PyTrilinos::Anasazi::EigensolverEpetra.


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

Generated on Thu Dec 17 11:00:20 2009 for PyTrilinos by  doxygen 1.5.9