I am trying to integrate a jump term over the interior boundary between two
3D elements. To do a sanity check on my code, I run the code snippet at the
end of this mail. It fails at the last assertion and I cannot figure out
why.
I am sorry if I am missing something simple.
Adam
/************* Example ********************/
// Assume that face has a non-null neighbor of elem.
// I am using MONOMIAL finite elements, TET4 elements and Fifth order
Gaussian quadrature.
// Both face_fe and neighbor_fe are initialized with an identical quadrature
rule.
const Elem* side_neighbor = elem->neighbor(face);
const unsigned int side_on_neighbor =
side_neighbor->which_neighbor_am_i(elem);
libmesh_assert(elem->neighbor(face) == side_neighbor);
libmesh_assert(side_neighbor->neighbor(side_on_neighbor) == elem);
face_fe->reinit(elem, face);
neighbor_fe->reinit(side_neighbor, side_on_neighbor);
const std::vector<Point>& f_point = face_fe->get_xyz();
const std::vector<Point>& n_point = neighbor_fe->get_xyz();
unsigned int f_size = f_point.size();
unsigned int n_size = n_point.size();
libmesh_assert(f_size == n_size);
for (unsigned int q_pt = 0; q_pt < f_size; ++q_pt)
{
const Point f = f_point[q_pt];
const Point n = n_point[q_pt];
libmesh_assert(f(0) == n(0) && f(1) == n(1) && f(2) == n(2)); //FAILS!
}
-------------------------------------------------------------------------
This SF.Net email is sponsored by the Moblin Your Move Developer's challenge
Build the coolest Linux based applications with Moblin SDK & win great prizes
Grand prize is a trip for two to an Open Source event anywhere in the world
http://moblin-contest.org/redirect.php?banner_id=100&url=/
_______________________________________________
Libmesh-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-users
