packages/Alea/mc_solvers/ChebyshevPolynomial.cc
//---------------------------------------------------------------------------//
// Compute smallest and largest magnitude eigenvalues
//---------------------------------------------------------------------------//
Compute smallest and largest magnitude eigenvalues
Compute smallest and largest magnitude eigenvalues
void ChebyshevPolynomial::computeEigenvalues()
{
    using std::chrono::high_resolution_clock;
    using std::chrono::duration_cast;
    using std::chrono::milliseconds;
    high_resolution_clock::time_point time_start = high_resolution_clock::now();
    if( b_verbosity >= HIGH )
        std::cout << "Computing matrix eigenvalues" << std::endl;
    bool zero_out = true;
    Teuchos::RCP<MV> x( new MV(b_A->getDomainMap(),1,zero_out) );
    typedef Anasazi::BasicEigenproblem<SCALAR,MV,OP> Eigenproblem;
    Teuchos::RCP<Eigenproblem> problem( new Eigenproblem(b_A,x) );
    problem->setNEV(1);
    problem->setProblem();
    Teuchos::ParameterList anasazi_pl;
    int global_length = x->getGlobalLength();
    int num_blocks = std::min(20,global_length-2);
    anasazi_pl.set<int>("Num Blocks",num_blocks);
    anasazi_pl.set<SCALAR>("Convergence Tolerance",1.0e-8);
    typedef Anasazi::BlockKrylovSchurSolMgr<SCALAR,MV,OP> KrylovSchur;
    anasazi_pl.set<std::string>("Which","LM");
    Teuchos::RCP<KrylovSchur> solver( new KrylovSchur(problem,anasazi_pl) );
    solver->solve();
    SCALAR realpart = solver->getRitzValues()[0].realpart;
    SCALAR imagpart = solver->getRitzValues()[0].imagpart;
    d_lambda_max = realpart;
    if( b_verbosity >= HIGH )
    {
        if( SCALAR_TRAITS::magnitude(imagpart) < 1.0e-8 )
        {
            std::cout << "Largest magnitude eigenvalue: "
                      << realpart << std::endl;
        }
        else
        {
            std::cout << "Largest magnitude eigenvalue: "
                      << realpart << " +/- " << imagpart << "i" << std::endl;
        }
    }
    anasazi_pl.set<std::string>("Which","SM");
    solver = Teuchos::rcp( new KrylovSchur(problem,anasazi_pl) );
    solver->solve();
    realpart = solver->getRitzValues()[0].realpart;
    imagpart = solver->getRitzValues()[0].imagpart;
    d_lambda_min = realpart;
    if( b_verbosity >= HIGH )
    {
        if( SCALAR_TRAITS::magnitude(imagpart) < 1.0e-8 )
        {
            std::cout << "Smallest magnitude eigenvalue: "
                      << realpart << std::endl;
        }
        else
        {
            std::cout << "Smallest magnitude eigenvalue: "
                      << realpart << " +/- " << imagpart << "i" << std::endl;
        }
    }
    high_resolution_clock::time_point time_end = high_resolution_clock::now();
    if( b_verbosity >= HIGH )
    {
        std::cout << "Eigenvalue calculation took "
            << duration_cast<milliseconds>(time_end-time_start).count()
            << " milliseconds" << std::endl;
    }
}
void ChebyshevPolynomial :: computeEigenvalues ( ) { using std :: chrono :: high_resolution_clock ; using std :: chrono :: duration_cast ; using std :: chrono :: milliseconds ; high_resolution_clock :: time_point time_start = high_resolution_clock :: now ( ) ; if ( b_verbosity >= HIGH ) std :: cout << " " << std :: endl ; bool zero_out = true ; Teuchos :: RCP < MV > x ( new MV ( b_A -> getDomainMap ( ) , 1 , zero_out ) ) ; typedef Anasazi :: BasicEigenproblem < SCALAR , MV , OP > Eigenproblem ; Teuchos :: RCP < Eigenproblem > problem ( new Eigenproblem ( b_A , x ) ) ; problem -> setNEV ( 1 ) ; problem -> setProblem ( ) ; Teuchos :: ParameterList anasazi_pl ; int global_length = x -> getGlobalLength ( ) ; int num_blocks = std :: min ( 20 , global_length - 2 ) ; anasazi_pl . set < int > ( " " , num_blocks ) ; anasazi_pl . set < SCALAR > ( " " , 1.0e-8 ) ; typedef Anasazi :: BlockKrylovSchurSolMgr < SCALAR , MV , OP > KrylovSchur ; anasazi_pl . set < std :: string > ( " " , " " ) ; Teuchos :: RCP < KrylovSchur > solver ( new KrylovSchur ( problem , anasazi_pl ) ) ; solver -> solve ( ) ; SCALAR realpart = solver -> getRitzValues ( ) [ 0 ] . realpart ; SCALAR imagpart = solver -> getRitzValues ( ) [ 0 ] . imagpart ; d_lambda_max = realpart ; if ( b_verbosity >= HIGH ) { if ( SCALAR_TRAITS :: magnitude ( imagpart ) < 1.0e-8 ) { std :: cout << " " << realpart << std :: endl ; } else { std :: cout << " " << realpart << " " << imagpart << " " << std :: endl ; } } anasazi_pl . set < std :: string > ( " " , " " ) ; solver = Teuchos :: rcp ( new KrylovSchur ( problem , anasazi_pl ) ) ; solver -> solve ( ) ; realpart = solver -> getRitzValues ( ) [ 0 ] . realpart ; imagpart = solver -> getRitzValues ( ) [ 0 ] . imagpart ; d_lambda_min = realpart ; if ( b_verbosity >= HIGH ) { if ( SCALAR_TRAITS :: magnitude ( imagpart ) < 1.0e-8 ) { std :: cout << " " << realpart << std :: endl ; } else { std :: cout << " " << realpart << " " << imagpart << " " << std :: endl ; } } high_resolution_clock :: time_point time_end = high_resolution_clock :: now ( ) ; if ( b_verbosity >= HIGH ) { std :: cout << " " << duration_cast < milliseconds > ( time_end - time_start ) . count ( ) << " " << std :: endl ; } }
Compute smallest and largest magnitude eigenvalues
Compute smallest and largest magnitude eigenvalues
// Create eigenvector // Anasazi eigenproblem // Set basic parameters // Use BlockKrylovSchur for now, can try out different solvers later // Compute largest magnitude // Compute largest magnitude