cmemory-corruptionpetsc

PETSc memory corruption on using DMDAVecGetArray


I'm trying to use PETSc's DMDA Vectors with 2 degrees of freedom and access them using struct, like in the manual. However, when I try to use DMDAVecGetArray even with one degree of freedom (like in example below) I get memory double free or corruption error. When I replace DMDAVecGetArray with VecGetArray everything works just fine.

What causes this error?


Compile MWE with

> gcc -Wall -Wextra -O0 -g $(pkg-config --cflags petsc mpi) -o mwe mwe.c $(pkg-config --libs petsc mpi)

and run with

> mpiexec --host localhost:$(nproc) -n $(nproc) mwe

MWE:

#include <petscdm.h>
#include <petscdmda.h>
#include <petscvec.h>

int
main(int argc, char **argv)
{
    PetscErrorCode ierr;
#define E(x)                                                                                                           \
    do {                                                                                                           \
        ierr = x;                                                                                              \
        CHKERRQ(ierr);                                                                                         \
    } while (0)

    PetscInt    N    = 10;
    Vec     fg, f;

    ierr = PetscInitialize(&argc, &argv, NULL, NULL);
    if (ierr) {
        return ierr;
    }

    DM domain;
    E(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_GHOSTED, N, 1, 1, NULL, &domain));
    E(DMSetUp(domain));

    E(DMCreateGlobalVector(domain, &fg));
    E(DMCreateLocalVector(domain, &f));
    E(VecZeroEntries(fg));

    E(VecAssemblyBegin(fg));
    E(VecAssemblyBegin(f));
    E(VecAssemblyEnd(fg));
    E(VecAssemblyEnd(f));

    PetscInt size;
    E(VecGetSize(f, &size));

    E(DMGlobalToLocal(domain, fg, INSERT_VALUES, f));

    PetscScalar *farr;
    E(VecGetArray(f, &farr)); // replace this
    // E(DMDAVecGetArrayWrite(domain, f, &farr)); // with this

    for (PetscInt i = 0; i < size; i++) {
        farr[i] += 1;
    }

    E(VecRestoreArray(f, &farr)); // replace this
    // E(DMDAVecRestoreArrayWrite(domain, f, &farr)); // with this
    E(DMLocalToGlobal(domain, f, INSERT_VALUES, fg));

    E(VecDestroy(&f));
    E(VecDestroy(&fg));
    E(DMDestroy(&domain));

    E(PetscFinalize());
    return ierr;
}

Debian bullseye; petsc-dev 3.14.4+dfsg1-1; libopenmpi-dev 4.1.0-7


Solution

  • Ok, found the problem. The key difference here is that when using VecGetArray underlying array is indexed using local indexes, and when using DMDAVecGetArray global indexes are used.

    So, in order to fix the code I have to get not only the size of the array, but also the first index of the local vector portion.

    Using function

    PetscInt start, size;
    E(DMDAGetCorners(domain, &start, NULL, NULL, &size, NULL, NULL));
    

    and then looping from start to start + size fixes the memory corruption error.