![]() |
Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include <linear_hex_map.hpp>
Public Member Functions | |
Linear_hex_map () | |
Linear_hex_map (const Linear_hex_map &) | |
template<typename Moab , typename Entity_handle , typename Points , typename Point > | |
std::pair< bool, CartVect > | evaluate_reverse (const double *verts, const double *eval_point, const double tol=1.e-6) const |
Private Member Functions | |
const CartVect & | reference_points (const std::size_t &i) const |
bool | is_contained (const CartVect ¶ms, const double tol) const |
bool | solve_inverse (const CartVect &point, CartVect ¶ms, const CartVect *verts, const double tol=1.e-6) const |
void | evaluate_forward (const CartVect &p, const CartVect *verts, CartVect &f) const |
double | integrate_scalar_field (const CartVect *points, const double *field_values) const |
template<typename Point , typename Points > | |
Matrix3 & | jacobian (const Point &p, const Points &points, Matrix3 &J) const |
Definition at line 16 of file linear_hex_map.hpp.
Definition at line 20 of file linear_hex_map.hpp.
{}
moab::element_utility::Linear_hex_map::Linear_hex_map | ( | const Linear_hex_map & | ) | [inline] |
Definition at line 22 of file linear_hex_map.hpp.
{}
void moab::element_utility::Linear_hex_map::evaluate_forward | ( | const CartVect & | p, |
const CartVect * | verts, | ||
CartVect & | f | ||
) | const [inline, private] |
Definition at line 107 of file linear_hex_map.hpp.
Referenced by moab::element_utility::Linear_hex_map< moab::Matrix3 >::solve_inverse().
{
typedef typename Points::value_type Vector;
f.set( 0.0, 0.0, 0.0 );
for( unsigned i = 0; i < 8; ++i )
{
const double N_i = ( 1 + p[0] * reference_points( i )[0] ) * ( 1 + p[1] * reference_points( i )[1] ) *
( 1 + p[2] * reference_points( i )[2] );
f += N_i * verts[i];
}
f *= 0.125;
return f;
}
std::pair< bool, CartVect > moab::element_utility::Linear_hex_map::evaluate_reverse | ( | const double * | verts, |
const double * | eval_point, | ||
const double | tol = 1.e-6 |
||
) | const [inline] |
Definition at line 27 of file linear_hex_map.hpp.
{
CartVect params( 3, 0.0 );
solve_inverse( eval_point, params, verts );
bool point_found = solve_inverse( eval_point, params, verts, tol ) && is_contained( params, tol );
return std::make_pair( point_found, params );
}
double moab::element_utility::Linear_hex_map::integrate_scalar_field | ( | const CartVect * | points, |
const double * | field_values | ||
) | const [inline, private] |
Definition at line 121 of file linear_hex_map.hpp.
{}
bool moab::element_utility::Linear_hex_map::is_contained | ( | const CartVect & | params, |
const double | tol | ||
) | const [inline, private] |
Definition at line 49 of file linear_hex_map.hpp.
Referenced by moab::element_utility::Linear_hex_map< moab::Matrix3 >::evaluate_reverse().
{
// just look at the box+tol here
return ( params[0] >= -1. - tol ) && ( params[0] <= 1. + tol ) && ( params[1] >= -1. - tol ) &&
( params[1] <= 1. + tol ) && ( params[2] >= -1. - tol ) && ( params[2] <= 1. + tol );
}
Matrix3& moab::element_utility::Linear_hex_map::jacobian | ( | const Point & | p, |
const Points & | points, | ||
Matrix3 & | J | ||
) | const [inline, private] |
Definition at line 124 of file linear_hex_map.hpp.
Referenced by moab::element_utility::Linear_hex_map< moab::Matrix3 >::solve_inverse().
{
}
const CartVect& moab::element_utility::Linear_hex_map::reference_points | ( | const std::size_t & | i | ) | const [inline, private] |
Definition at line 41 of file linear_hex_map.hpp.
Referenced by moab::element_utility::Linear_hex_map< moab::Matrix3 >::evaluate_forward().
{
const CartVect rpts[8] = { CartVect( -1, -1, -1 ), CartVect( 1, -1, -1 ), CartVect( 1, 1, -1 ),
CartVect( -1, 1, -1 ), CartVect( -1, -1, 1 ), CartVect( 1, -1, 1 ),
CartVect( 1, 1, 1 ), CartVect( -1, 1, 1 ) };
return rpts[i];
}
bool moab::element_utility::Linear_hex_map::solve_inverse | ( | const CartVect & | point, |
CartVect & | params, | ||
const CartVect * | verts, | ||
const double | tol = 1.e-6 |
||
) | const [inline, private] |
Definition at line 56 of file linear_hex_map.hpp.
Referenced by moab::element_utility::Linear_hex_map< moab::Matrix3 >::evaluate_reverse().
{
const double error_tol_sqr = tol * tol;
CartVect delta( 0.0, 0.0, 0.0 );
params = delta;
evaluate_forward( params, verts, delta );
delta -= point;
std::size_t num_iterations = 0;
#ifdef LINEAR_HEX_DEBUG
std::stringstream ss;
ss << "CartVect: ";
ss << point[0] << ", " << point[1] << ", " << point[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( delta.length_squared() > error_tol_sqr )
{
#ifdef LINEAR_HEX_DEBUG
ss << "Iter #: " << num_iterations << " Err: " << delta.length() << " Iterate: ";
ss << params[0] << ", " << params[1] << ", " << params[2] << std::endl;
#endif
if( ++num_iterations >= 5 )
{
return false;
}
Matrix3 J;
jacobian( params, verts, J );
double det = moab::Matrix3::determinant3( J );
if( fabs( det ) < 1.e-10 )
{
#ifdef LINEAR_HEX_DEBUG
std::cerr << ss.str();
#endif
#ifndef LINEAR_HEX_DEBUG
std::cerr << x[0] << ", " << x[1] << ", " << x[2] << std::endl;
#endif
std::cerr << "inverse solve failure: det: " << det << std::endl;
exit( -1 );
}
params -= moab::Matrix3::inverse( J, 1.0 / det ) * delta;
evaluate_forward( params, points, delta );
delta -= x;
}
return true;
}