idaholab / moose

Multiphysics Object Oriented Simulation Environment
https://www.mooseframework.org
GNU Lesser General Public License v2.1
1.71k stars 1.04k forks source link

Three-dimensional contact robustness issues #16840

Open recuero opened 3 years ago

recuero commented 3 years ago

Bug Description

A three-dimensional problem with second-order tets results in segmentation fault when the bodies come into contact, as reported in #16773.

Steps to Reproduce

The problem built in #16773 reproduces this problem.

Impact

Reduced robustness of three-dimensional Tet10 (possibly Tet4) mechanical contact.

recuero commented 3 years ago

For the model reported with second-order tet meshes, when the bodies come into contact, there is a failing libmesh check. This makes simulation stop also in opt mode. Backtrace below.

    frame #1: 0x00000001088e4dc6 libmesh_dbg.0.dylib`libMesh::libmesh_terminate_handler() at libmesh.C:321:5
    frame #2: 0x00007fff6bea3887 libc++abi.dylib`std::__terminate(void (*)()) + 8
    frame #3: 0x00007fff6bea61a2 libc++abi.dylib`__cxxabiv1::failed_throw(__cxxabiv1::__cxa_exception*) + 27
    frame #4: 0x00007fff6bea6169 libc++abi.dylib`__cxa_throw + 113
    frame #5: 0x00000001002441a0 libcontact-dbg.0.dylib`libMesh::DofObject::dof_number(this=0x0000000156daf358, s=1, var=11, comp=0) const at dof_object.h:993:3
    frame #6: 0x000000010027c2d8 libcontact-dbg.0.dylib`MechanicalContactConstraint::nodalArea(this=0x000000018931e618, pinfo=0x00000001a81e7150) at MechanicalContactConstraint.C:1724:27
    frame #7: 0x000000010027825f libcontact-dbg.0.dylib`MechanicalContactConstraint::getPenalty(this=0x000000018931e618, pinfo=0x00000001a81e7150) at MechanicalContactConstraint.C:1743:16
    frame #8: 0x0000000100279174 libcontact-dbg.0.dylib`MechanicalContactConstraint::computeContactForce(this=0x000000018931e618, pinfo=0x00000001a81e7150, update_contact_set=true) at MechanicalContactConstraint.C:556:24
    frame #9: 0x0000000100282c14 libcontact-dbg.0.dylib`MechanicalContactConstraint::shouldApply(this=0x000000018931e618) at MechanicalContactConstraint.C:499:9
    frame #10: 0x0000000102f8cee7 libmoose-dbg.0.dylib`NonlinearSystemBase::constraintResiduals(this=0x0000000156016e18, residual=0x0000000167c4d1e0, displaced=true) at NonlinearSystemBase.C:1179:24
    frame #11: 0x0000000102f86844 libmoose-dbg.0.dylib`NonlinearSystemBase::computeResidualInternal(this=0x0000000156016e18, tags=size=3) at NonlinearSystemBase.C:1586:9
    frame #12: 0x0000000102f85022 libmoose-dbg.0.dylib`NonlinearSystemBase::computeResidualTags(this=0x0000000156016e18, tags=size=3) at NonlinearSystemBase.C:732:5
    frame #13: 0x00000001027d01ca libmoose-dbg.0.dylib`FEProblemBase::computeResidualTags(this=0x0000000156014818, tags=size=3) at FEProblemBase.C:5524:8
    frame #14: 0x00000001027cf776 libmoose-dbg.0.dylib`FEProblemBase::computeResidualInternal(this=0x0000000156014818, soln=0x0000000114c5b360, residual=0x00007ffeefbfe3d0, tags=size=3) at FEProblemBase.C:5362:5
    frame #15: 0x00000001027cf59c libmoose-dbg.0.dylib`FEProblemBase::computeResidual(this=0x0000000156014818, soln=0x0000000114c5b360, residual=0x00007ffeefbfe3d0) at FEProblemBase.C:5317:3
    frame #16: 0x00000001027f24d6 libmoose-dbg.0.dylib`FEProblemBase::computeResidualSys(this=0x0000000156014818, (null)=0x0000000114c5ad40, soln=0x0000000114c5b360, residual=0x00007ffeefbfe3d0) at FEProblemBase.C:5292:3
    frame #17: 0x0000000102fa1807 libmoose-dbg.0.dylib`ComputeResidualFunctor::residual(this=0x0000000156018310, soln=0x0000000114c5b360, residual=0x00007ffeefbfe3d0, sys=0x0000000114c5ad40) at ComputeResidualFunctor.C:23:15
    frame #18: 0x00000001096591ce libmesh_dbg.0.dylib`::libmesh_petsc_snes_residual(snes=0x0000000218813620, x=0x0000000218809020, r=0x00000001b000be20, ctx=0x0000000114c5b9f0) at petsc_nonlinear_solver.C:163:35
    frame #19: 0x000000010e6fcfae libpetsc.3.14.dylib`SNESComputeFunction + 334
    frame #20: 0x000000010e711bd4 libpetsc.3.14.dylib`SNESSolve_NEWTONLS + 404
    frame #21: 0x000000010e70253d libpetsc.3.14.dylib`SNESSolve + 1837
    frame #22: 0x000000010965e051 libmesh_dbg.0.dylib`libMesh::PetscNonlinearSolver<double>::solve(this=0x0000000114c5b9f0, pre_in=0x0000000156ba1380, x_in=0x0000000114c5b240, r_in=0x0000000114c5b810, (null)=0.001, (null)=25) at petsc_nonlinear_solver.C:951:10
    frame #23: 0x00000001097034c8 libmesh_dbg.0.dylib`libMesh::NonlinearImplicitSystem::solve(this=0x0000000114c5ad40) at nonlinear_implicit_system.C:178:27
    frame #24: 0x0000000103291d6d libmoose-dbg.0.dylib`TimeIntegrator::solve(this=0x0000000167ce6398) at TimeIntegrator.C:66:16
    frame #25: 0x0000000102fa28f8 libmoose-dbg.0.dylib`NonlinearSystem::solve(this=0x0000000156016e18) at NonlinearSystem.C:200:23
    frame #26: 0x00000001027fa7c0 libmoose-dbg.0.dylib`FEProblemBase::solve(this=0x0000000156014818) at FEProblemBase.C:5042:10
    frame #27: 0x00000001022da134 libmoose-dbg.0.dylib`FEProblemSolve::solve(this=0x0000000188010980) at FEProblemSolve.C:197:12
    frame #28: 0x00000001022d5496 libmoose-dbg.0.dylib`PicardSolve::solveStep(this=0x0000000188010a68, begin_norm_old=1.7976931348623157E+308, begin_norm=0x0000000167c60840, end_norm_old=1.7976931348623157E+308, end_norm=0x0000000167c9c160, relax=false, relaxed_dofs=size=0) at PicardSolve.C:445:22
    frame #29: 0x00000001022d38c7 libmoose-dbg.0.dylib`PicardSolve::solve(this=0x0000000188010a68) at PicardSolve.C:236:28
    frame #30: 0x0000000102d89b9c libmoose-dbg.0.dylib`TimeStepper::step(this=0x0000000167d3c8a8) at TimeStepper.C:163:43
    frame #31: 0x00000001022d7ea0 libmoose-dbg.0.dylib`Transient::takeStep(this=0x0000000188010818, input_dt=-1) at Transient.C:431:20
    frame #32: 0x00000001022e08fb libmoose-dbg.0.dylib`Transient::execute(this=0x0000000188010818) at Transient.C:314:5
    frame #33: 0x00000001035d8cd8 libmoose-dbg.0.dylib`MooseApp::executeExecutioner(this=0x0000000157008e18) at MooseApp.C:1023:19
    frame #34: 0x00000001035d40bf libmoose-dbg.0.dylib`MooseApp::run(this=0x0000000157008e18) at MooseApp.C:1150:5
    frame #35: 0x0000000100002f96 contact-dbg`main(argc=3, argv=0x00007ffeefbffa18) at main.C:36:8

Retrieving the dof number appears to be ultimately causing the issue:

dof_id_type dof = node->dof_number(_aux_system.number(), _nodal_area_var->number(), 0);

sapitts commented 3 years ago

I ran into this problem with a Quad8 mesh as well, so it looks like the issue is not just confined to tet elements. Stack trace is quite similar:

Assertion `comp < this->n_comp(s,var)' failed.
comp = 0
this->n_comp(s,var) = 0

Stack frames: 36
0: 0   libmesh_dbg.0.dylib                 0x0000000115aa9e0c libMesh::print_trace(std::__1::basic_ostream<char, std::__1::char_traits<char> >&) + 460
1: 1   libmesh_dbg.0.dylib                 0x0000000115a9ec6a libMesh::MacroFunctions::report_error(char const*, int, char const*, char const*) + 218
2: 2   libcontact-dbg.0.dylib              0x000000010861ef43 libMesh::DofObject::dof_number(unsigned int, unsigned int, unsigned int) const + 1587
3: 3   libcontact-dbg.0.dylib              0x00000001086586d8 MechanicalContactConstraint::nodalArea(PenetrationInfo&) + 104
4: 4   libcontact-dbg.0.dylib              0x0000000108657db4 MechanicalContactConstraint::computeContactForce(PenetrationInfo*, bool) + 12132
5: 5   libcontact-dbg.0.dylib              0x000000010865f145 MechanicalContactConstraint::shouldApply() + 245
6: 6   libmoose-dbg.0.dylib                0x000000010feb5c60 NonlinearSystemBase::constraintResiduals(libMesh::NumericVector<double>&, bool) + 1408
7: 7   libmoose-dbg.0.dylib                0x000000010feafa17 NonlinearSystemBase::computeResidualInternal(std::__1::set<unsigned int, std::__1::less<unsigned int>, std::__1::allocator<unsigned int> > const&) + 4055
8: 8   libmoose-dbg.0.dylib                0x000000010feae3d2 NonlinearSystemBase::computeResidualTags(std::__1::set<unsigned int, std::__1::less<unsigned int>, std::__1::allocator<unsigned int> > const&) + 1282
9: 9   libmoose-dbg.0.dylib                0x000000010f6b61da FEProblemBase::computeResidualTags(std::__1::set<unsigned int, std::__1::less<unsigned int>, std::__1::allocator<unsigned int> > const&) + 2346
10: 10  libmoose-dbg.0.dylib                0x000000010f6b5786 FEProblemBase::computeResidualInternal(libMesh::NumericVector<double> const&, libMesh::NumericVector<double>&, std::__1::set<unsigned int, std::__1::less<unsigned int>, std::__1::allocator<unsigned int> > const&) + 230
11: 11  libmoose-dbg.0.dylib                0x000000010f6b558c FEProblemBase::computeResidual(libMesh::NumericVector<double> const&, libMesh::NumericVector<double>&) + 252
12: 12  libmoose-dbg.0.dylib                0x000000010f6d8f36 FEProblemBase::computeResidualSys(libMesh::NonlinearImplicitSystem&, libMesh::NumericVector<double> const&, libMesh::NumericVector<double>&) + 102
13: 13  libmoose-dbg.0.dylib                0x000000010feca947 ComputeResidualFunctor::residual(libMesh::NumericVector<double> const&, libMesh::NumericVector<double>&, libMesh::NonlinearImplicitSystem&) + 87
14: 14  libmesh_dbg.0.dylib                 0x000000011678bcae libmesh_petsc_snes_residual + 1678
15: 15  libpetsc.3.14.2.dylib               0x00000001192c555e SNESComputeFunction + 334
16: 16  libpetsc.3.14.2.dylib               0x000000011929eae4 SNESLineSearchApply_Basic + 388
17: 17  libpetsc.3.14.2.dylib               0x0000000119297607 SNESLineSearchApply + 263
18: 18  libpetsc.3.14.2.dylib               0x00000001192e899b SNESSolve_NEWTONLS + 2395
19: 19  libpetsc.3.14.2.dylib               0x00000001192cabcd SNESSolve + 1837
20: 20  libmesh_dbg.0.dylib                 0x00000001167908f6 libMesh::PetscNonlinearSolver<double>::solve(libMesh::SparseMatrix<double>&, libMesh::NumericVector<double>&, libMesh::NumericVector<double>&, double, unsigned int) + 3926
21: 21  libmesh_dbg.0.dylib                 0x000000011682d9d4 libMesh::NonlinearImplicitSystem::solve() + 596
22: 22  libmoose-dbg.0.dylib                0x000000011020962d TimeIntegrator::solve() + 45
23: 23  libmoose-dbg.0.dylib                0x000000010fecba38 NonlinearSystem::solve() + 1016
24: 24  libmoose-dbg.0.dylib                0x000000010f6e1080 FEProblemBase::solve() + 384
25: 25  libmoose-dbg.0.dylib                0x000000010f18f9b4 FEProblemSolve::solve() + 36
26: 26  libmoose-dbg.0.dylib                0x000000010f18aac6 PicardSolve::solveStep(double, double&, double, double&, bool, std::__1::set<unsigned long long, std::__1::less<unsigned long long>, std::__1::allocator<unsigned long long> > const&) + 1622
27: 27  libmoose-dbg.0.dylib                0x000000010f188c60 PicardSolve::solve() + 3728
28: 28  libmoose-dbg.0.dylib                0x000000010fcae48c TimeStepper::step() + 44
29: 29  libmoose-dbg.0.dylib                0x000000010f18d610 Transient::takeStep(double) + 288
30: 30  libmoose-dbg.0.dylib                0x000000010f19617b Transient::execute() + 171
31: 31  libmoose-dbg.0.dylib                0x000000011055c368 MooseApp::executeExecutioner() + 232
32: 32  libmoose-dbg.0.dylib                0x00000001105576cf MooseApp::run() + 1647
33: 33  freya-dbg                           0x000000010721dfb5 main + 149
34: 34  libdyld.dylib                       0x00007fff68ec2cc9 start + 1
35: 35  ???                                 0x0000000000000003 0x0 + 3
[0] /Users/pittsa/miniconda3/envs/moose/libmesh/include/libmesh/dof_object.h, line 993, compiled Mar 22 2021 at 17:14:44