On Mon, Jul 18, 2016 at 9:42 AM, Xujun Zhao <xzha...@gmail.com> wrote:
> Guys, > > I met this error when I try to determine if a point is inside an element > using the function constains_point(). I tracked the code using the > dbg version and found that this happened when computing the > FE::inverse_map(), in which it calls TypeTensor::inverse() and leads to > Jacobian = 0. Here are the error messages: > > Assertion `my_det != static_cast<T>(0.)' failed. > > my_det = 0 > > static_cast<T>(0.) = 0 > > > Stack frames: 12 > > 0: 0 libmesh_dbg.0.dylib 0x000000010b72db6c > libMesh::print_trace(std::ostream&) + 472 > > 1: 1 libmesh_dbg.0.dylib 0x000000010b729122 > libMesh::MacroFunctions::report_error(char const*, int, char const*, char > const*) + 387 > > 2: 2 libmesh_dbg.0.dylib 0x000000010bde82f2 > libMesh::TypeTensor<double>::inverse() const + 666 > > 3: 3 libmesh_dbg.0.dylib 0x000000010b93c35a > libMesh::FE<3u, (libMesh::FEFamily)0>::inverse_map(libMesh::Elem const*, > libMesh::Point const&, double, bool) + 1250 > > 4: 4 libmesh_dbg.0.dylib 0x000000010b8a101d > libMesh::FEInterface::inverse_map(unsigned int, libMesh::FEType const&, > libMesh::Elem const*, libMesh::Point const&, double, bool) + 3725 > > 5: 5 libmesh_dbg.0.dylib 0x000000010bae6f18 > libMesh::Elem::point_test(libMesh::Point const&, double, double) const + > 1744 > > 6: 6 libmesh_dbg.0.dylib 0x000000010bae680a > libMesh::Elem::contains_point(libMesh::Point const&, double) const + 250 > > 7: 7 example-dbg 0x000000010b2e7fd2 > libMesh::PointMesh<3u>::build_elem_neighbor_list() + 1158 > > 8: 8 example-dbg 0x000000010b2e8aae > libMesh::PointMesh<3u>::reinit() + 534 > > 9: 9 example-dbg 0x000000010b275e07 > test_cell(libMesh::Parallel::Communicator const&) + 4146 > > > I am curious why this happens. It seems to me that J can be zero when the > volume of an element is zero. But I let elem->volume() to compute the > volume value of each element, none of those are zeros! > > Also, I am using the second order Hex27 element, can this be caused by the > unconvergence of newton iteration? Thanks a lot. > Yeah... the Newton iteration inside FE::inverse_map() appears to have failed, I suppose this could happen if the Newton method evaluates the Jacobian at a point outside the element, but I'm a little surprised the determinant was *exactly* 0.0. BTW, a non-zero volume() doesn't guarantee a pointwise non-singular Jacobian (especially for a complicated element like the Hex27). It just depends on where one is evaluating it, i.e. quadrature points vs. arbitrary points within the element. It would be helpful if you could post the node positions of the element in question (in the libmesh node ordering for Hex27). That way we could take a closer look. -- John ------------------------------------------------------------------------------ What NetFlow Analyzer can do for you? Monitors network bandwidth and traffic patterns at an interface-level. Reveals which users, apps, and protocols are consuming the most bandwidth. Provides multi-vendor support for NetFlow, J-Flow, sFlow and other flows. Make informed decisions using capacity planning reports.http://sdm.link/zohodev2dev _______________________________________________ Libmesh-users mailing list Libmesh-users@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/libmesh-users