[petsc-users] Getting a vector from a DM to output VTK

[Matthew Knepley]

On Thu, Dec 22, 2022 at 12:41 AM Nicholas Arnold-Medabalimi <
narnoldm at umich.edu> wrote:

> Hi Petsc Users
>
> I've been having trouble consistently getting a vector generated from a DM
> to output to VTK correctly. I've used ex1.c (which works properly)to try
> and figure it out, but I'm still having some issues. I must be missing
> something small that isn't correctly associating the section with the DM.
>
>     DMPlexGetChart(dm, &p0, &p1);
>     PetscSection section_full;
>     PetscSectionCreate(PETSC_COMM_WORLD, &section_full);
>     PetscSectionSetNumFields(section_full, 1);
>     PetscSectionSetChart(section_full, p0, p1);
>     PetscSectionSetFieldName(section_full, 0, "state");
>
>     for (int i = c0; i < c1; i++)
>     {
>         PetscSectionSetDof(section_full, i, 1);
>         PetscSectionSetFieldDof(section_full, i, 0, 1);
>     }
>     PetscSectionSetUp(section_full);
>     DMSetNumFields(dm, 1);
>     DMSetLocalSection(dm, section_full);
>     DMCreateGlobalVector(dm, &state_full);
>
>     int o0, o1;
>     VecGetOwnershipRange(state_full, &o0, &o1);
>     PetscScalar *state_full_array;
>     VecGetArray(state_full, &state_full_array);
>
>     for (int i = 0; i < (c1 - c0); i++)
>     {
>         int offset;
>         PetscSectionGetOffset(section_full, i, &offset);
>         state_full_array[offset] = 101325 + i;
>     }
>
>     VecRestoreArray(state_full, &state_full_array);
>
>
>     PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
>     PetscViewerSetType(viewer, PETSCVIEWERVTK);
>     PetscViewerFileSetMode(viewer, FILE_MODE_WRITE);
>     PetscViewerFileSetName(viewer, "mesh.vtu");
>     VecView(state_full, viewer);
>
> If I run this mesh.vtu isn't generated at all. If I instead do a DMView
> passing the DM it will just output the mesh correctly.
>
> Any assistance would be greatly appreciated.
>

DMCreateGlobalVector() dispatches to DMCreateGlobalVector_Plex(), which
resets the view method to VecView_Plex(), which should dispatch to
VecView_Plex_Local_VTK(). You can verify this in the debugger, or send us
code we can run to verify it.

  Thanks,

    Matt


> Sincerely
> Nicholas
>
>
> --
> Nicholas Arnold-Medabalimi
>
> Ph.D. Candidate
> Computational Aeroscience Lab
> University of Michigan
>


-- 
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/20221223/9571aff3/attachment.html>