MOAB: Mesh Oriented datABase  (version 5.3.1)
test_Matrix3.cpp File Reference
#include <vector>
#include <cassert>
#include "moab/Matrix3.hpp"
#include "moab/Core.hpp"
#include "moab/CartVect.hpp"
#include "TestUtil.hpp"
+ Include dependency graph for test_Matrix3.cpp:

Go to the source code of this file.

Defines

#define ACOS(x)   acos( std::min( std::max( x, -1.0 ), 1.0 ) )
#define CHECK_EIGVECREAL_EQUAL(EXP, ACT, EPS)   check_equal_eigvect( ( EXP ), ( ACT ), ( EPS ), #EXP, #ACT, __LINE__, __FILE__ )

Functions

double find_angle (const moab::CartVect &A, const moab::CartVect &B)
void check_equal_eigvect (const moab::CartVect &A, const moab::CartVect &B, double eps, const char *sA, const char *sB, int line, const char *file)
void test_EigenDecomp ()
int main ()

Define Documentation

#define ACOS (   x)    acos( std::min( std::max( x, -1.0 ), 1.0 ) )

Definition at line 10 of file test_Matrix3.cpp.

Referenced by find_angle().

#define CHECK_EIGVECREAL_EQUAL (   EXP,
  ACT,
  EPS 
)    check_equal_eigvect( ( EXP ), ( ACT ), ( EPS ), #EXP, #ACT, __LINE__, __FILE__ )

Definition at line 22 of file test_Matrix3.cpp.

Referenced by test_EigenDecomp().


Function Documentation

void check_equal_eigvect ( const moab::CartVect A,
const moab::CartVect B,
double  eps,
const char *  sA,
const char *  sB,
int  line,
const char *  file 
)

Definition at line 24 of file test_Matrix3.cpp.

References moab::angle(), check_equal(), eps, find_angle(), flag_error(), and moab::CartVect::length().

{
    check_equal( A.length(), B.length(), eps, sA, sB, line, file );

    double angle = find_angle( A, B );

    if( ( fabs( A[0] - B[0] ) <= eps || fabs( A[0] + B[0] ) <= eps ) &&
        ( fabs( A[1] - B[1] ) <= eps || fabs( A[1] + B[1] ) <= eps ) &&
        ( fabs( A[2] - B[2] ) <= eps || fabs( A[2] + B[2] ) <= eps ) &&
        ( angle <= eps || fabs( angle - 180.0 ) <= eps ) )
        return;

    std::cout << "Equality Test Failed: " << sA << " == " << sB << std::endl;
    std::cout << "  at line " << line << " of '" << file << "'" << std::endl;

    std::cout << "  Expected: ";
    std::cout << A << std::endl;

    std::cout << "  Actual:   ";
    std::cout << B << std::endl;

    std::cout << "Angle between vectors := " << angle << std::endl;

    flag_error();
}
double find_angle ( const moab::CartVect A,
const moab::CartVect B 
)

Definition at line 12 of file test_Matrix3.cpp.

References ACOS, moab::dot(), and moab::CartVect::length().

Referenced by check_equal_eigvect().

{
    const double lA = A.length(), lB = B.length();
    assert( lA > 0.0 );
    assert( lB > 0.0 );
    const double dPI = 3.14159265;
    const double dot = ( A[0] * B[0] + A[1] * B[1] + A[2] * B[2] );
    return ACOS( dot / ( lA * lB ) ) * 180.0 / dPI;
}
int main ( )

Definition at line 53 of file test_Matrix3.cpp.

References PointInVol::result, RUN_TEST, and test_EigenDecomp().

{

    int result = 0;

    result += RUN_TEST( test_EigenDecomp );

    return result;
}
void test_EigenDecomp ( )

Definition at line 65 of file test_Matrix3.cpp.

References CHECK_EIGVECREAL_EQUAL, CHECK_ERR, CHECK_REAL_EQUAL, moab::Matrix3::col(), moab::Matrix3::eigen_decomposition(), ErrorCode, moab::CartVect::length(), moab::Matrix::matrix_vector(), and moab::CartVect::normalize().

Referenced by main().

{
    // Create a matrix
    moab::Matrix3 mat;

    mat( 0 ) = 2;
    mat( 1 ) = -1;
    mat( 2 ) = 0;
    mat( 3 ) = -1;
    mat( 4 ) = 2;
    mat( 5 ) = -1;
    mat( 6 ) = 0;
    mat( 7 ) = -1;
    mat( 8 ) = 2;

    // now do the Eigen Decomposition of this Matrix

    CartVect lamda;
    Matrix3 vectors;
    moab::ErrorCode rval = mat.eigen_decomposition( lamda, vectors );CHECK_ERR( rval );
    for( int i = 0; i < 3; ++i )
        vectors.col( i ).normalize();

    // Hardcoded check values for the results
    double lamda_check[3];
    lamda_check[0] = 0.585786;
    lamda_check[1] = 2.0;
    lamda_check[2] = 3.41421;

    moab::CartVect vec0_check( 0.5, 0.707107, 0.5 );
    moab::CartVect vec1_check( 0.707107, 3.37748e-17, -0.707107 );
    moab::CartVect vec2_check( 0.5, -0.707107, 0.5 );

    // now verfy that the returns Eigenvalues and Eigenvectors are correct (within some tolerance)
    double tol = 1e-04;

    // use a deprecated constructor
    Matrix3 mat3( CartVect( 2, -1, 0 ), CartVect( -1, 2, -1 ), CartVect( 0, -1, 2 ), true );
    CHECK_REAL_EQUAL( mat( 1 ), mat3( 1 ), tol );

    vec0_check.normalize();
    vec1_check.normalize();
    vec2_check.normalize();

    // check that the correct Eigenvalues are returned correctly (in order)
    CHECK_REAL_EQUAL( lamda[0], lamda_check[0], tol );
    CHECK_REAL_EQUAL( lamda[1], lamda_check[1], tol );
    CHECK_REAL_EQUAL( lamda[2], lamda_check[2], tol );

    // check the Eigenvector values (order should correspond to the Eigenvalues)
    // first vector
    CHECK_EIGVECREAL_EQUAL( vectors.col( 0 ), vec0_check, tol );

    // sceond vector
    CHECK_EIGVECREAL_EQUAL( vectors.col( 1 ), vec1_check, tol );

    // third vector
    CHECK_EIGVECREAL_EQUAL( vectors.col( 2 ), vec2_check, tol );

    // another check to ensure the result is valid (AM-kM = 0)
    for( unsigned i = 0; i < 3; ++i )
    {
        moab::CartVect v = moab::Matrix::matrix_vector( mat, vectors.col( i ) ) - lamda[i] * vectors.col( i );
        CHECK_REAL_EQUAL( v.length(), 0, tol );
    }

    // for a real, symmetric matrix the Eigenvectors should be orthogonal
    CHECK_REAL_EQUAL( vectors.col( 0 ) % vectors.col( 1 ), 0, tol );
    CHECK_REAL_EQUAL( vectors.col( 0 ) % vectors.col( 2 ), 0, tol );

    return;
}
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines