<div>Hello,</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">I ran into a seemingly non-consistent behavior of the function <span></span><pre><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/DMPlex/DMPlexComputeCellGeometryFEM/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_rJABo0_$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">DMPlexComputeCellGeometryFEM</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">(</span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/DM/DM/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_kyYWLFi$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">DM</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">dm</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">,</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/Sys/PetscInt/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_jAVrwAd$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscInt</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">cell</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">,</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/DT/PetscQuadrature/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_jVjnRvX$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscQuadrature</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">quad</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">,</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/Sys/PetscReal/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_heQHCRx$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscReal</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">*</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">v</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">,</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/Sys/PetscReal/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_heQHCRx$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscReal</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">J</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">[],</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/Sys/PetscReal/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_heQHCRx$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscReal</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">invJ</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">[],</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span><a href="https://urldefense.us/v3/__https://petsc.org/release/manualpages/Sys/PetscReal/__;!!G_uCfscf7eWS!YEv7_aaCdN0cxV8sUhRBIACYZGNNcwQusNjEpfVep6yOAxRYE9f8kw4XiVENrjEvQEVwiheebMTQVhPw4eenYApT_heQHCRx$"><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">PetscReal</span></a></span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);"> </span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">*</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">detJ</span><span style="font-family: Menlo, Consolas, Courier New, monospace; color: rgb(0, 0, 0);">)</span></pre>Calling it in a loop happens to crash, only some times, when the input <span style="font-family: Menlo, Consolas, Courier New, monospace;">npoints</span> is not equal to one.</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Attached is a mesh file (<span style="font-family: Menlo, Consolas, Courier New, monospace;">msh</span> format), consisting of four quadrilaterals. Small examples in both C and Fortran also attached; they read the mesh, create a Gauss tensor quadrature rule and then loop through all cells in order to get the mapped quadrature points.<div><div><div><div><div><span><br></span></div><div><span>When the quadrature created has a single point per dimension (<span style="font-family: Menlo, Consolas, Courier New, monospace;">npoints = 1</span>), there seems to be no issue. The program runs and the output is the expected. </span></div><div><span><br></span></div><div><span>On the C side, For npoints = 2 the program hangs, usually no errors; sometimes</span></div><div><span><br></span></div><div><span> <span style="font-family: Menlo, Consolas, Courier New, monospace;">*** stack smashing detected ***: terminated</span>,</span></div><div><span><span style="font-family: Menlo, Consolas, Courier New, monospace;">[:1345101] *** Process received signal ***</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">[:1345101] Signal: Aborted (6)</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">[:1345101] Signal code:  (-6)</span></div><span></span><br></span></div><div><span>or</span></div><div><span><br></span></div><div><span> <span style="font-family: Menlo, Consolas, Courier New, monospace;">LeakSanitizer: detected memory leaks</span> (when compiled with</span><span style="font-family: Menlo, Consolas, Courier New, monospace;"> -fsanitize=addres).</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><br></span></div><div><span style="font-family: Arial, sans-serif;">On the Fortran side, similar behavior for </span><span style="font-family: Menlo, Consolas, Courier New, monospace;">npoints = 1 or 2</span><span style="font-family: Arial, sans-serif;">. </span><span style="font-family: Arial, sans-serif;">Either no error at all (no error meaning the one that only mentions "</span><span style="font-family: Menlo, Consolas, Courier New, monospace;">probably memory access out of bounds</span><span style="font-family: Arial, sans-serif;">"), or sometimes "</span><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span style="font-family: Menlo, Consolas, Courier New, monospace;">corrupted size vs. prev_size</span></span><span style="font-family: Arial, sans-serif;">", followed by a stack trace. </span></div><div><span style="font-family: Arial, sans-serif;"><br></span></div><div><span style="font-family: Arial, sans-serif;">Question: Do the arguments </span><span style="font-family: Menlo, Consolas, Courier New, monospace;">J, invJ </span><span style="font-family: Arial, sans-serif;">need to be allocated (e.g. </span><span style="font-family: Menlo, Consolas, Courier New, monospace;">allocate(J(4))</span><span style="font-family: Arial, sans-serif;">) before calling </span><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span style="font-family: Menlo, Consolas, Courier New, monospace;">DMPlexComputeCellGeometryFEM</span></span><span style="font-family: Arial, sans-serif;">? Allocation seems to be the cause for triggering an error.</span><span style="font-family: Menlo, Consolas, Courier New, monospace;"><br></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><br></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">---<br></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><br></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">NOTE:</span><span style="font-family: Arial, sans-serif;"> The interface to </span><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span><span style="font-family: Menlo, Consolas, Courier New, monospace;">DMPlexComputeCellGeometryFEM</span><span style="font-family: Arial, sans-serif;">, under </span><span style="font-family: Menlo, Consolas, Courier New, monospace;">$PETSC_ARCH$/ftn/dm/petscdmplex.h90</span><span style="font-family: Arial, sans-serif;">, read:</span></span></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span><span style="font-family: Arial, sans-serif;"><br></span></span></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span><span style="font-family: Arial, sans-serif;"><span style="font-family: Menlo, Consolas, Courier New, monospace;">  interface DMPlexComputeCellGeometryFEM</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  subroutine DMPlexComputeCellGeometryFEM(a,b,c,d,e,f,g, z)</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  import tDM,tPetscQuadrature</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  DM :: a</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscInt :: b</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscQuadrature :: c</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscReal :: d          ####   should it be d(*) ?<br></span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscReal :: e(*)</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscReal :: f(*)</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscReal :: g</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  PetscErrorCode z</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  end subroutine</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">  end interface</span></div><span></span><br></span></span><span style="font-family:Arial, sans-serif"></span></span></div></div></div></div></div></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Should argument <span style="font-family: Menlo, Consolas, Courier New, monospace;">d</span> (referring to the array of the mapped integration points <span style="font-family: Menlo, Consolas, Courier New, monospace;">v</span>) be defined differently? <br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---<br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Arial, sans-serif;">NOTE: Tested with PETSc 3.23.1 (from the tarball) and the latest git.</span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Arial, sans-serif;">Configured with</span><span style="font-family: Menlo, Consolas, Courier New, monospace;"> ./configure --with-fortran-bindings=1</span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;">Using gcc/gfortran 11.4.0, OpenMPI 4.1.2<br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;"><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;">Thanks,<br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;">Noam<br></span></div>