Power Method Examples
[Assorted Thyra Operator/Vector Example Code]

Here we show some simple examples of using code in the Thyra package with an example linear ANA algorithm for the power method for estimating the dominate eigen value of a matrix. More...

Modules
group	Templated Serial Implementation of the Power Method
	Here is an example program that shows the use of the example serial templated matrix class ExampleTridiagSerialLinearOp with the example linear ANA implementation sillyPowerMethod().
Functions
template<class Scalar>
bool	sillyPowerMethod (const Thyra::LinearOpBase< Scalar > &A, const int maxNumIters, const typename Teuchos::ScalarTraits< Scalar >::magnitudeType tolerance, Scalar *lambda, std::ostream *out=NULL)
	Silly little example power method abstract numerical algorithm (ANA).

---

Detailed Description

Here we show some simple examples of using code in the Thyra package with an example linear ANA algorithm for the power method for estimating the dominate eigen value of a matrix.

These example programs are meant to mimic the power method example program shown in the Epetra documentation.

The power method ANA is implemented in the function sillyPowerMethod() and its implementation is shown below:

template<class Scalar>
bool sillyPowerMethod(
  const Thyra::LinearOpBase<Scalar>                              &A
  ,const int                                                     maxNumIters
  ,const typename Teuchos::ScalarTraits<Scalar>::magnitudeType   tolerance
  ,Scalar                                                        *lambda
  ,std::ostream                                                  *out          = NULL
  )
{
  // Create some typedefs and some other stuff to make the code cleaner
  typedef Teuchos::ScalarTraits<Scalar> ST; typedef typename ST::magnitudeType ScalarMag;
  const Scalar one = ST::one(); using Thyra::NOTRANS;
  typedef Teuchos::RefCountPtr<const Thyra::VectorSpaceBase<Scalar> > VectorSpacePtr;
  typedef Teuchos::RefCountPtr<Thyra::VectorBase<Scalar> > VectorPtr;
  // Initialize
  if(out) *out << "\nStarting power method (target tolerrance = "<<tolerance<<") ...\n\n";
  VectorPtr q = createMember(A.domain()), z = createMember(A.range()), r = createMember(A.range());
  Thyra::seed_randomize<Scalar>(0);
  Thyra::randomize( Scalar(-one), Scalar(+one), &*z );
  // Perform iterations
  for( int iter = 0; iter < maxNumIters; ++iter ) {
    const ScalarMag z_nrm = norm(*z);       // Compute natural norm of z
    V_StV( &*q, Scalar(one/z_nrm), *z );    // q = (1/||z}*z 
    apply( A, NOTRANS , *q, &*z );          // z = A*q
    *lambda = scalarProd(*q,*z);            // lambda = <q,z>    : Approximate maximum absolute eigenvalue
    if( iter%(maxNumIters/10) == 0 || iter+1 == maxNumIters ) {
      V_StVpV(&*r,Scalar(-*lambda),*q,*z);  // r = -lambda*q + z : Compute residual of eigenvalue equation
      const ScalarMag r_nrm = norm(*r);     // Compute natural norm of r
      if(out) *out << "Iter = " << iter << ", lambda = " << (*lambda) << ", ||A*q-lambda*q|| = " << r_nrm << std::endl;
      if( r_nrm < tolerance ) return true;  // Success!
    }
  }
  if(out) *out << "\nMaximum number of iterations exceeded with ||-lambda*q + z|| > tolerence = " << tolerance << "\n";
  return false; // Failure
} // end sillyPowerMethod

The above templated function sillyPowerMethod() is used in the following various example implementations which use several different scalar types:

1. Templated Serial Implementation of the Power Method

---

Function Documentation

template<class Scalar> bool sillyPowerMethod (  const Thyra::LinearOpBase< Scalar > &  A, const int  maxNumIters, const typename Teuchos::ScalarTraits< Scalar >::magnitudeType  tolerance, Scalar *  lambda, std::ostream *  out = NULL )

Silly little example power method abstract numerical algorithm (ANA). This little function is just a silly little ANA that implements the power method for estimating the dominate eigenvalue for a matrix using the foundational Thyra operator/vector interfaces. This function is small and is meant to be looked at so study its implementation by clicking on the below link to its definition. Examples: sillyPowerMethod.hpp, and sillyPowerMethod_serial.cpp. Definition at line 47 of file sillyPowerMethod.hpp.

---