Actual source code: ex169.c

  1: static char help[] = "Test memory leak when duplicating a redundant matrix.\n\n";

  3: /*
  4:   Include "petscmat.h" so that we can use matrices.
  5:   automatically includes:
  6:      petscsys.h    - base PETSc routines   petscvec.h    - vectors
  7:      petscmat.h    - matrices
  8:      petscis.h     - index sets            petscviewer.h - viewers
  9: */
 10: #include <petscmat.h>

 12: int main(int argc, char **args)
 13: {
 14:   Mat          A, Ar, C;
 15:   PetscViewer  fd;                       /* viewer */
 16:   char         file[PETSC_MAX_PATH_LEN]; /* input file name */
 17:   PetscInt     ns = 2;
 18:   PetscMPIInt  size;
 19:   PetscSubcomm subc;
 20:   PetscBool    flg;

 22:   PetscFunctionBeginUser;
 23:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 24:   /*
 25:      Determine files from which we read the two linear systems
 26:      (matrix and right-hand-side vector).
 27:   */
 28:   PetscCall(PetscOptionsGetString(NULL, NULL, "-f0", file, sizeof(file), &flg));
 29:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate binary file with the -f0 option");
 30:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
 31:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 32:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Reading matrix with %d processors\n", size));
 33:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 34:   PetscCall(MatLoad(A, fd));
 35:   PetscCall(PetscViewerDestroy(&fd));
 36:   /*
 37:      Determines amount of subcomunicators
 38:   */
 39:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nsub", &ns, NULL));
 40:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Splitting in %" PetscInt_FMT " subcommunicators\n", ns));
 41:   PetscCall(PetscSubcommCreate(PetscObjectComm((PetscObject)A), &subc));
 42:   PetscCall(PetscSubcommSetNumber(subc, ns));
 43:   PetscCall(PetscSubcommSetType(subc, PETSC_SUBCOMM_CONTIGUOUS));
 44:   PetscCall(PetscSubcommSetFromOptions(subc));
 45:   PetscCall(MatCreateRedundantMatrix(A, 0, PetscSubcommChild(subc), MAT_INITIAL_MATRIX, &Ar));
 46:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Copying matrix\n"));
 47:   PetscCall(MatDuplicate(Ar, MAT_COPY_VALUES, &C));
 48:   PetscCall(MatAXPY(Ar, 0.1, C, DIFFERENT_NONZERO_PATTERN));
 49:   PetscCall(PetscSubcommDestroy(&subc));

 51:   /*
 52:      Free memory
 53:   */
 54:   PetscCall(MatDestroy(&A));
 55:   PetscCall(MatDestroy(&Ar));
 56:   PetscCall(MatDestroy(&C));
 57:   PetscCall(PetscFinalize());
 58:   return 0;
 59: }

 61: /*TEST

 63:    test:
 64:       nsize: 4
 65:       requires: !complex double !defined(PETSC_USE_64BIT_INDICES)
 66:       args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/ns-real-int32-float64 -malloc_dump

 68: TEST*/