Content

Eigenvalues (Symmetric Matrix)

In this example we compute the eigenvalues \(w \in \mathbb{C}\) and eigenvectors \(V = (v_1, \dots, v_n)\) of a symmetric Matrix \(A \in \mathbb{R}\).

We use lapack::ev which is the FLENS port of LAPACK's dsyev.

Example Code

#include <iostream>
#include <flens/flens.cxx>

using namespace std;
using namespace flens;

typedef double   T;

int
main()
{
    typedef SyMatrix<FullStorage<T, ColMajor> >   SYMatrix;
    typedef DenseVector<Array<T> >                Vector;

    const int n = 5;

    SYMatrix  A(n, Lower);
    Vector    w(n);

    A.diag( 0) =  2, -5,  60, -6;
    A.diag(-1) =  3,  0, -68;
    A.diag(-2) = -1,  2,  2;
    A.diag(-3) =  0, -6;
    A.diag(-4) =  2;

    cerr << "A = " << A << endl;

    Vector   work;

    // int     optSize = ev_wsq(true, A);
    // Vector  work(optSize);

    lapack::ev(true, A, w, work);

    cerr << "Eigenvalues: w = " << w << endl;
    cerr << "Eigenvectors: V = " << A.general() << endl;
}

Comments on Example Code

Define some convenient typedefs for the matrix/vector types of our system of linear equations.

    typedef SyMatrix<FullStorage<T, ColMajor> >   SYMatrix;
    typedef DenseVector<Array<T> >                Vector;

Vector for workspace. If this vector has zero length then lapack::ev will do a worksize query and also resize work.

    Vector   work;

You also could do a worksize query manually

    // int     optSize = ev_wsq(true, A);
    // Vector  work(optSize);

Call lapack::ev to compute eigenvalues \(w\) and eigenvectors \(V\). Note that the (full-)storage reserved by the symmetric matrix \(A\) (i.e. the lower and upper part) gets overwritten with the eigenvectors.

    lapack::ev(true, A, w, work);

We use the general view (A.general()) for accessing the eigenvectors stored in the underlying full storage of \(A\):

    cerr << "Eigenvalues: w = " << w << endl;
    cerr << "Eigenvectors: V = " << A.general() << endl;
}

Compile and Run

Currently no generic implementation of dsyev exists. So we have to link against a native LAPACK implementation:

$shell> cd flens/examples                                                      
$shell> g++ -std=c++11 -Wall -I../.. -DWITH_VECLIB -DUSE_CXXLAPACK -framework vecLib -o lapack-syev lapack-syev.cc                    
$shell> ./lapack-syev                                                          
A = 
            2             3            -1             0             2 
            3            -5             0             2            -6 
           -1             0             6            -6             2 
            0             2            -6             0             8 
            2            -6             2             8            -6 
Eigenvalues: w = 
     -16.5476       -3.49574        2.81132        3.81611        10.4159 
Eigenvectors: V = 
    -0.163881      0.286144      0.901433     -0.237839      0.148713 
     0.484838     -0.738709      0.207981     -0.418743     0.0252737 
    -0.192772     -0.317136      0.307277       0.44326     -0.755887 
    -0.465057     -0.495194     0.0821765      0.412236      0.601509 
     0.696158      0.163209      0.207341      0.633773      0.209923