[petsc-users] DM/DMPlex Issue from DM vector access

[Matthew Knepley]

On Mon, Dec 6, 2021 at 6:06 PM Ferrand, Jesus A. <FERRANJ2 at my.erau.edu>
wrote:

> Dear PETSc Team:
>
> I am a new PETSc user who is working on an FEA code and ran into an issue
> pertaining to DMPlex.
>
> I have a gmsh mesh file that I import using "DMPlexCreateGmshFromFile()."
> I then fetch the XYZ coordinates of the nodes from this mesh using
> "DMGetCoordinatesLocal()." Deeper into the code, I have a call to
> "DMPlexGetTransitiveClosure()" inside a loop that scans the cells (I think
> you guys refer to it as "depth-3")
>

It is more robust to use height 0 for cells.


> to reference nodes. Here's the catch: If I attempt to reference the array
> from the vector (Vec) in the call to "DMGetCoordinatesLocal()" by using
> something like "VecGetArray()," or "VecGetArrayRead()," the call to
> "DMPlexGetTransitiveClosure()" errors-out with a segmentation fault. I need
> access to that vector's stored XYZ-data because I'm using my own finite
> element scripts. I have no clue as to why this is happening. Maybe it is
> a newbie mistake and I am forgetting to restore some memory?
>

1) The First DMCreate() is unnecessary and leaks right now

2) PetscSectionDestroy(&s) after DMSetLocalSection()

3) Once you have set the section, you can use DMCreateGlobalVector() and
DMCreateMatrix() instead of doing it by hand

4) DMPlexGetTransitiveClosure() has somewhat complicated memory management:

  a) The safest thing to do is initialize closure to NULL, and call
Get/RestoreClosure() each iteration

  b) If you want to manage the memory yourself, then allocate closure[]
coming in and every iteration reset size_closure to the array size on input

  c) You can combine these, so that you set closure to NULL on input, so
that GetClosure() allocates the array for you. Then make the rest of your
       calls without changing size_closure and closure. Then after the loop
ends call RestoreClosure() and it will deallocate.

Note that this looks like elasticity. SNES ex17 is my idea of doing
elasticity :)

  Thanks,

     Matt


> Code, error message, and gmsh files are attached.
> Your help is much appreciated.
>
> Machine Type: HP Laptop
> C-compiler: Gnu C
> OS: Ubuntu 20.04
> PETSc version: 3.16.0
> MPI Implementation: MPICH
>
> Sincerely:
>
> *J.A. Ferrand*
>
> Embry-Riddle Aeronautical University - Daytona Beach FL
>
> M.Sc. Aerospace Engineering | May 2022
>
> B.Sc. Aerospace Engineering
>
> B.Sc. Computational Mathematics
>
>
>
> Sigma Gamma Tau
>
> Tau Beta Pi
>
> Honors Program
>
>
>
> *Phone:* (386)-843-1829
>
> *Email(s):* ferranj2 at my.erau.edu
>
>     jesus.ferrand at gmail.com
>


-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20211206/255d403a/attachment.html>