dune-grid
2.2.0
|
00001 #ifndef DUNE_ALBERTA_TRANSFORMATION_HH 00002 #define DUNE_ALBERTA_TRANSFORMATION_HH 00003 00004 #include <dune/common/fvector.hh> 00005 00006 #include <dune/grid/albertagrid/misc.hh> 00007 00008 #if HAVE_ALBERTA 00009 00010 namespace Dune 00011 { 00012 00013 class AlbertaTransformation 00014 { 00015 typedef Alberta::GlobalSpace GlobalSpace; 00016 00017 public: 00018 typedef Alberta::Real ctype; 00019 00020 static const int dimension = Alberta::dimWorld; 00021 00022 typedef FieldVector< ctype, dimension > WorldVector; 00023 00024 explicit 00025 AlbertaTransformation ( const Alberta::AffineTransformation *trafo = NULL ) 00026 : matrix_( (trafo != NULL ? trafo->M : GlobalSpace::identityMatrix()) ), 00027 shift_( (trafo != NULL ? trafo->t : GlobalSpace::nullVector()) ) 00028 {} 00029 00030 AlbertaTransformation ( const GlobalSpace::Matrix &matrix, 00031 const GlobalSpace::Vector &shift ) 00032 : matrix_( matrix ), 00033 shift_( shift ) 00034 {} 00035 00036 WorldVector evaluate ( const WorldVector &x ) const 00037 { 00038 WorldVector y; 00039 for( int i = 0; i < dimension; ++i ) 00040 { 00041 const GlobalSpace::Vector &row = matrix_[ i ]; 00042 y[ i ] = shift_[ i ]; 00043 for( int j = 0; j < dimension; ++j ) 00044 y[ i ] += row[ j ] * x[ j ]; 00045 } 00046 return y; 00047 } 00048 00049 WorldVector evaluateInverse ( const WorldVector &y ) const 00050 { 00051 // Note: ALBERTA requires the matrix to be orthogonal 00052 WorldVector x( ctype( 0 ) ); 00053 for( int i = 0; i < dimension; ++i ) 00054 { 00055 const GlobalSpace::Vector &row = matrix_[ i ]; 00056 const ctype v = y[ i ] - shift_[ i ]; 00057 for( int j = 0; j < dimension; ++j ) 00058 x[ j ] += row[ j ] * v; 00059 } 00060 return x; 00061 } 00062 00063 private: 00064 const GlobalSpace::Matrix &matrix_; 00065 const GlobalSpace::Vector &shift_; 00066 }; 00067 00068 } 00069 00070 #endif // #if HAVE_ALBERTA 00071 00072 #endif // #ifndef DUNE_ALBERTA_TRANSFORMATION_HH