![]() |
Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include <quadratic_hex_map.hpp>
Public Member Functions | |
Quadratic_hex_map () | |
Quadratic_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) const |
Static Public Member Functions | |
typedef _Matrix | Matrix |
Private Types | |
typedef Quadratic_hex_map< Matrix > | Self |
Private Member Functions | |
double | reference_points (const std::size_t &i, const std::size_t &j) const |
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) const |
template<typename Point , typename Points > | |
Point & | evaluate (const Point &p, const Points &points, Point &f) const |
template<typename Point , typename Field > | |
double | evaluate_scalar_field (const Point &p, const Field &field) const |
template<typename Field , typename Points > | |
double | integrate_scalar_field (const Points &p, const Field &field_values) const |
template<typename Point , typename Points > | |
Matrix & | jacobian (const Point &p, const Points &, Matrix &J) const |
Definition at line 51 of file quadratic_hex_map.hpp.
typedef _Matrix moab::element_utility::Quadratic_hex_map< _Matrix >::Matrix [static] |
Definition at line 54 of file quadratic_hex_map.hpp.
typedef Quadratic_hex_map< Matrix > moab::element_utility::Quadratic_hex_map< _Matrix >::Self [private] |
Definition at line 57 of file quadratic_hex_map.hpp.
moab::element_utility::Quadratic_hex_map< _Matrix >::Quadratic_hex_map | ( | ) | [inline] |
Definition at line 61 of file quadratic_hex_map.hpp.
{}
moab::element_utility::Quadratic_hex_map< _Matrix >::Quadratic_hex_map | ( | const Self & | f | ) | [inline] |
Definition at line 63 of file quadratic_hex_map.hpp.
{}
Point& moab::element_utility::Quadratic_hex_map< _Matrix >::evaluate | ( | const Point & | p, |
const Points & | points, | ||
Point & | f | ||
) | const [inline, private] |
Definition at line 169 of file quadratic_hex_map.hpp.
Referenced by moab::element_utility::Quadratic_hex_map< moab::Matrix3 >::solve_inverse().
{
typedef typename Points::value_type Vector;
Vector result;
for( int i = 0; i < 3; ++i )
{
result[i] = 0;
}
for( unsigned i = 0; i < 27; ++i )
{
const double sh = SH( reference_points( i, 0 ), p[0] ) * SH( reference_points( i, 1 ), p[1] ) *
SH( reference_points( i, 2 ), p[2] );
result += sh * points[i];
}
for( int i = 0; i < 3; ++i )
{
f[i] = result[i];
}
return f;
}
double moab::element_utility::Quadratic_hex_map< _Matrix >::evaluate_scalar_field | ( | const Point & | p, |
const Field & | field | ||
) | const [inline, private] |
Definition at line 190 of file quadratic_hex_map.hpp.
{
double x = 0.0;
for( int i = 0; i < 27; i++ )
{
const double sh = SH( reference_points( i, 0 ), p[0] ) * SH( reference_points( i, 1 ), p[1] ) *
SH( reference_points( i, 2 ), p[2] );
x += sh * field[i];
}
return x;
}
double moab::element_utility::Quadratic_hex_map< _Matrix >::integrate_scalar_field | ( | const Points & | p, |
const Field & | field_values | ||
) | const [inline, private] |
Definition at line 202 of file quadratic_hex_map.hpp.
{
// TODO: gaussian integration , probably 2x2x2
return 0.;
}
bool moab::element_utility::Quadratic_hex_map< _Matrix >::is_contained | ( | const Point & | p, |
const double | tol | ||
) | const [inline, private] |
Definition at line 110 of file quadratic_hex_map.hpp.
Referenced by moab::element_utility::Quadratic_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 );
}
Matrix& moab::element_utility::Quadratic_hex_map< _Matrix >::jacobian | ( | const Point & | p, |
const Points & | , | ||
Matrix & | J | ||
) | const [inline, private] |
Definition at line 209 of file quadratic_hex_map.hpp.
Referenced by moab::element_utility::Quadratic_hex_map< moab::Matrix3 >::solve_inverse().
{
J = Matrix( 0.0 );
for( int i = 0; i < 27; i++ )
{
const double sh[3] = { SH( reference_points( i, 0 ), p[0] ), SH( reference_points( i, 1 ), p[1] ),
SH( reference_points( i, 2 ), p[2] ) };
const double dsh[3] = { DSH( reference_points( i, 0 ), p[0] ), DSH( reference_points( i, 1 ), p[1] ),
DSH( reference_points( i, 2 ), p[2] ) };
for( int j = 0; j < 3; j++ )
{
// dxj/dr first column
J( j, 0 ) += dsh[0] * sh[1] * sh[2] * reference_points( i, j );
J( j, 1 ) += sh[0] * dsh[1] * sh[2] * reference_points( i, j ); // dxj/ds
J( j, 2 ) += sh[0] * sh[1] * dsh[2] * reference_points( i, j ); // dxj/dt
}
}
return J;
}
std::pair< bool, Point > moab::element_utility::Quadratic_hex_map< _Matrix >::operator() | ( | const Moab & | , |
const Entity_handle & | , | ||
const Points & | v, | ||
const Point & | p, | ||
const double | tol = 1.e-6 |
||
) | const [inline] |
Definition at line 68 of file quadratic_hex_map.hpp.
{
Point result( 3, 0.0 );
bool point_found = solve_inverse( p, result, v, tol ) && is_contained( result, tol );
return std::make_pair( point_found, result );
}
double moab::element_utility::Quadratic_hex_map< _Matrix >::reference_points | ( | const std::size_t & | i, |
const std::size_t & | j | ||
) | const [inline, private] |
Definition at line 83 of file quadratic_hex_map.hpp.
Referenced by moab::element_utility::Quadratic_hex_map< moab::Matrix3 >::evaluate(), moab::element_utility::Quadratic_hex_map< moab::Matrix3 >::evaluate_scalar_field(), and moab::element_utility::Quadratic_hex_map< moab::Matrix3 >::jacobian().
{
const double rpts[27][3] = { { -1, -1, -1 }, { 1, -1, -1 }, { 1, 1, -1 }, // reference_points nodes: 0-7
{ -1, 1, -1 }, // mid-edge nodes: 8-19
{ -1, -1, 1 }, // center-face nodes 20-25 center node 26
{ 1, -1, 1 }, //
{ 1, 1, 1 }, { -1, 1, 1 }, // 4 ----- 19 ----- 7
{ 0, -1, -1 }, // . | . |
{ 1, 0, -1 }, // 16 25 18 |
{ 0, 1, -1 }, // . | . |
{ -1, 0, -1 }, // 5 ----- 17 ----- 6 |
{ -1, -1, 0 }, // | 12 | 23 15
{ 1, -1, 0 }, // | | |
{ 1, 1, 0 }, // | 20 | 26 | 22 |
{ -1, 1, 0 }, // | | |
{ 0, -1, 1 }, // 13 21 | 14 |
{ 1, 0, 1 }, // | 0 ----- 11 ----- 3
{ 0, 1, 1 }, // | . | .
{ -1, 0, 1 }, // | 8 24 | 10
{ 0, -1, 0 }, // | . | .
{ 1, 0, 0 }, // 1 ----- 9 ----- 2
{ 0, 1, 0 }, //
{ -1, 0, 0 }, { 0, 0, -1 }, { 0, 0, 1 }, { 0, 0, 0 } };
return rpts[i][j];
}
bool moab::element_utility::Quadratic_hex_map< _Matrix >::solve_inverse | ( | const Point & | x, |
Point & | xi, | ||
const Points & | points, | ||
const double | tol = 1.e-6 |
||
) | const [inline, private] |
Definition at line 118 of file quadratic_hex_map.hpp.
Referenced by moab::element_utility::Quadratic_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 QUADRATIC_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 QUADRATIC_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 QUADRATIC_HEX_DEBUG
std::cerr << ss.str();
#endif
#ifndef QUADRATIC_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;
}