MOAB: Mesh Oriented datABase  (version 5.4.1)
moab::element_utility::Spectral_hex_map< _Matrix > Class Template Reference

#include <spectral_hex_map.hpp>

+ Inheritance diagram for moab::element_utility::Spectral_hex_map< _Matrix >:
+ Collaboration diagram for moab::element_utility::Spectral_hex_map< _Matrix >:

Public Member Functions

 Spectral_hex_map ()
 Spectral_hex_map (int order)
 Spectral_hex_map (const Self &f)
template<typename Moab , typename Entity_handle , typename Points , typename Point >
std::pair< bool, Point > operator() (const Moab &, const Entity_handle &, const Points &v, const Point &p, const double tol=1.e-6)

Static Public Member Functions

typedef _Matrix Matrix

Private Types

typedef Spectral_hex_map< MatrixSelf

Private Member Functions

void initialize_spectral_hex (int order)
void free_data ()
void set_gl_points (double *x, double *y, double *z)
template<typename Point >
bool is_contained (const Point &p, const double tol) const
template<typename Point , typename Points >
bool solve_inverse (const Point &x, Point &xi, const Points &points, const double tol=1.e-6)
template<typename Point , typename Points >
Point & evaluate (const Point &p, const Points &, Point &f)
template<typename Point , typename Field >
double evaluate_scalar_field (const Point &p, const Field &field) const
template<typename Points , typename Field >
double integrate_scalar_field (const Points &p, const Field &field) const
template<typename Point , typename Points >
Matrixjacobian (const Point &, const Points &, Matrix &J)

Private Attributes

bool _init
int _n
real_z [3]
lagrange_data _ld [3]
opt_data_3 _data
real_odwork
real_xyz [3]

Detailed Description

template<typename _Matrix>
class moab::element_utility::Spectral_hex_map< _Matrix >

Definition at line 22 of file spectral_hex_map.hpp.


Member Typedef Documentation

template<typename _Matrix>
typedef _Matrix moab::element_utility::Spectral_hex_map< _Matrix >::Matrix [static]

Definition at line 25 of file spectral_hex_map.hpp.

template<typename _Matrix>
typedef Spectral_hex_map< Matrix > moab::element_utility::Spectral_hex_map< _Matrix >::Self [private]

Definition at line 28 of file spectral_hex_map.hpp.


Constructor & Destructor Documentation

template<typename _Matrix>
moab::element_utility::Spectral_hex_map< _Matrix >::Spectral_hex_map ( ) [inline]

Definition at line 32 of file spectral_hex_map.hpp.

{};
template<typename _Matrix>
moab::element_utility::Spectral_hex_map< _Matrix >::Spectral_hex_map ( int  order) [inline]

Definition at line 33 of file spectral_hex_map.hpp.

        {
            initialize_spectral_hex( order );
        }
template<typename _Matrix>
moab::element_utility::Spectral_hex_map< _Matrix >::Spectral_hex_map ( const Self f) [inline]

Definition at line 38 of file spectral_hex_map.hpp.

{}

Member Function Documentation

template<typename _Matrix>
template<typename Point , typename Points >
Point& moab::element_utility::Spectral_hex_map< _Matrix >::evaluate ( const Point &  p,
const Points ,
Point &  f 
) [inline, private]

Definition at line 160 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::solve_inverse().

        {
            for( int d = 0; d < 3; ++d )
            {
                lagrange_0( &_ld[d], p[0] );
            }
            for( int d = 0; d < 3; ++d )
            {
                f[d] = tensor_i3( _ld[0].J, _ld[0].n, _ld[1].J, _ld[1].n, _ld[2].J, _ld[2].n, _xyz[d], _odwork );
            }
            return f;
        }
template<typename _Matrix>
template<typename Point , typename Field >
double moab::element_utility::Spectral_hex_map< _Matrix >::evaluate_scalar_field ( const Point &  p,
const Field &  field 
) const [inline, private]

Definition at line 174 of file spectral_hex_map.hpp.

        {
            int d;
            for( d = 0; d < 3; d++ )
            {
                lagrange_0( &_ld[d], p[d] );
            }
            return tensor_i3( _ld[0].J, _ld[0].n, _ld[1].J, _ld[1].n, _ld[2].J, _ld[2].n, field, _odwork );
        }
template<typename _Matrix>
void moab::element_utility::Spectral_hex_map< _Matrix >::free_data ( ) [inline, private]

Definition at line 63 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::initialize_spectral_hex().

        {
            for( int d = 0; d < 3; d++ )
            {
                free( _z[d] );
                lagrange_free( &_ld[d] );
            }
            opt_free_3( &_data );
            free( _odwork );
        }
template<typename _Matrix>
void moab::element_utility::Spectral_hex_map< _Matrix >::initialize_spectral_hex ( int  order) [inline, private]

Definition at line 41 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::Spectral_hex_map().

        {
            if( _init && _n == order )
            {
                return;
            }
            if( _init && _n != order )
            {
                free_data();
            }
            _init = true;
            _n    = order;
            for( int d = 0; d < 3; d++ )
            {
                lobatto_nodes( _z[d], _n );
                lagrange_setup( &_ld[d], _z[d], _n );
            }
            opt_alloc_3( &_data, _ld );
            std::size_t nf = _n * _n, ne = _n, nw = 2 * _n * _n + 3 * _n;
            _odwork = tmalloc( real, 6 * nf + 9 * ne + nw );
        }
template<typename _Matrix>
template<typename Points , typename Field >
double moab::element_utility::Spectral_hex_map< _Matrix >::integrate_scalar_field ( const Points p,
const Field &  field 
) const [inline, private]

Definition at line 184 of file spectral_hex_map.hpp.

        {
            // set the position of GL points
            // set the positions of GL nodes, before evaluations
            _data.elx[0] = _xyz[0];
            _data.elx[1] = _xyz[1];
            _data.elx[2] = _xyz[2];
            double xi[3];
            // triple loop; the most inner loop is in r direction, then s, then t
            double integral = 0.;
            // double volume = 0;
            int index = 0;  // used fr the inner loop
            for( int k = 0; k < _n; k++ )
            {
                xi[2] = _ld[2].z[k];
                // double wk= _ld[2].w[k];
                for( int j = 0; j < _n; j++ )
                {
                    xi[1] = _ld[1].z[j];
                    // double wj= _ld[1].w[j];
                    for( int i = 0; i < _n; i++ )
                    {
                        xi[0] = _ld[0].z[i];
                        // double wi= _ld[0].w[i];
                        opt_vol_set_intp_3( (opt_data_3*)&_data, xi );  // cast away const-ness
                        double wk = _ld[2].J[k];
                        double wj = _ld[1].J[j];
                        double wi = _ld[0].J[i];
                        Matrix3 J( 0. );
                        for( int n = 0; n < 8; n++ )
                            J( n / 3, n % 3 ) = _data.jac[n];
                        double bm = wk * wj * wi * J.determinant();
                        integral += bm * field[index++];
                        // volume +=bm;
                    }
                }
            }
            // std::cout << "volume: " << volume << "\n";
            return integral;
        }
template<typename _Matrix>
template<typename Point >
bool moab::element_utility::Spectral_hex_map< _Matrix >::is_contained ( const Point &  p,
const double  tol 
) const [inline, private]

Definition at line 101 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::operator()().

        {
            // just look at the box+tol here
            return ( p[0] >= -1. - tol ) && ( p[0] <= 1. + tol ) && ( p[1] >= -1. - tol ) && ( p[1] <= 1. + tol ) &&
                   ( p[2] >= -1. - tol ) && ( p[2] <= 1. + tol );
        }
template<typename _Matrix>
template<typename Point , typename Points >
Matrix& moab::element_utility::Spectral_hex_map< _Matrix >::jacobian ( const Point &  ,
const Points ,
Matrix J 
) [inline, private]

Definition at line 226 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::solve_inverse().

        {
            real x[3];
            for( int i = 0; i < 3; ++i )
            {
                _data.elx[i] = _xyz[i];
            }
            opt_vol_set_intp_3( &_data, x );
            for( int i = 0; i < 9; ++i )
            {
                J( i % 3, i / 3 ) = _data.jac[i];
            }
            return J;
        }
template<typename _Matrix>
template<typename Moab , typename Entity_handle , typename Points , typename Point >
std::pair< bool, Point > moab::element_utility::Spectral_hex_map< _Matrix >::operator() ( const Moab &  ,
const Entity_handle &  ,
const Points v,
const Point &  p,
const double  tol = 1.e-6 
) [inline]

Definition at line 77 of file spectral_hex_map.hpp.

        {
            Point result( 3, 0.0 );
            /*
            moab.tag_get_by_ptr(_xm1Tag, &eh, 1,(const void **) &_xyz[ 0] );
            moab.tag_get_by_ptr(_ym1Tag, &eh, 1,(const void **) &_xyz[ 1] );
            moab.tag_get_by_ptr(_zm1Tag, &eh, 1,(const void **) &_xyz[ 2] );
            */
            bool point_found = solve_inverse( p, result, v, tol ) && is_contained( result, tol );
            return std::make_pair( point_found, result );
        }
template<typename _Matrix>
void moab::element_utility::Spectral_hex_map< _Matrix >::set_gl_points ( double *  x,
double *  y,
double *  z 
) [inline, private]

Definition at line 94 of file spectral_hex_map.hpp.

        {
            _xyz[0] = x;
            _xyz[1] = y;
            _xyz[2] = z;
        }
template<typename _Matrix>
template<typename Point , typename Points >
bool moab::element_utility::Spectral_hex_map< _Matrix >::solve_inverse ( const Point &  x,
Point &  xi,
const Points points,
const double  tol = 1.e-6 
) [inline, private]

Definition at line 109 of file spectral_hex_map.hpp.

Referenced by moab::element_utility::Spectral_hex_map< moab::Matrix3 >::operator()().

        {
            const double error_tol_sqr = tol * tol;
            Point delta( 3, 0.0 );
            xi = delta;
            evaluate( xi, points, delta );
            vec_subtract( delta, x );
            std::size_t num_iterations = 0;
#ifdef SPECTRAL_HEX_DEBUG
            std::stringstream ss;
            ss << "Point: ";
            ss << x[0] << ", " << x[1] << ", " << x[2] << std::endl;
            ss << "Hex: ";
            for( int i = 0; i < 8; ++i )
            {
                ss << points[i][0] << ", " << points[i][1] << ", " << points[i][2] << std::endl;
            }
            ss << std::endl;
#endif
            while( normsq( delta ) > error_tol_sqr )
            {
#ifdef SPECTRAL_HEX_DEBUG
                ss << "Iter #: " << num_iterations << " Err: " << sqrt( normsq( delta ) ) << " Iterate: ";
                ss << xi[0] << ", " << xi[1] << ", " << xi[2] << std::endl;
#endif
                if( ++num_iterations >= 5 )
                {
                    return false;
                }
                Matrix J;
                jacobian( xi, points, J );
                double det = moab::Matrix::determinant3( J );
                if( fabs( det ) < 1.e-10 )
                {
#ifdef SPECTRAL_HEX_DEBUG
                    std::cerr << ss.str();
#endif
#ifndef SPECTRAL_HEX_DEBUG
                    std::cerr << x[0] << ", " << x[1] << ", " << x[2] << std::endl;
#endif
                    std::cerr << "inverse solve failure: det: " << det << std::endl;
                    exit( -1 );
                }
                vec_subtract( xi, moab::Matrix::inverse( J, 1.0 / det ) * delta );
                evaluate( xi, points, delta );
                vec_subtract( delta, x );
            }
            return true;
        }

Member Data Documentation

template<typename _Matrix>
bool moab::element_utility::Spectral_hex_map< _Matrix >::_init [private]

List of all members.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines