Actual source code: ex169.c


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

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

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

 23:   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:   PetscOptionsGetString(NULL,NULL,"-f0",file,sizeof(file),&flg);
 30:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
 31:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 32:   PetscPrintf(PETSC_COMM_WORLD,"Reading matrix with %d processors\n",size);
 33:   MatCreate(PETSC_COMM_WORLD,&A);
 34:   MatLoad(A,fd);
 35:   PetscViewerDestroy(&fd);
 36:   /*
 37:      Determines amount of subcomunicators
 38:   */
 39:   PetscOptionsGetInt(NULL,NULL,"-nsub",&ns,NULL);
 40:   PetscPrintf(PETSC_COMM_WORLD,"Splitting in %" PetscInt_FMT " subcommunicators\n",ns);
 41:   PetscSubcommCreate(PetscObjectComm((PetscObject)A),&subc);
 42:   PetscSubcommSetNumber(subc,ns);
 43:   PetscSubcommSetType(subc,PETSC_SUBCOMM_CONTIGUOUS);
 44:   PetscSubcommSetFromOptions(subc);
 45:   MatCreateRedundantMatrix(A,0,PetscSubcommChild(subc),MAT_INITIAL_MATRIX,&Ar);
 46:   PetscPrintf(PETSC_COMM_WORLD,"Copying matrix\n");
 47:   MatDuplicate(Ar,MAT_COPY_VALUES,&C);
 48:   MatAXPY(Ar,0.1,C,DIFFERENT_NONZERO_PATTERN);
 49:   PetscSubcommDestroy(&subc);

 51:   /*
 52:      Free memory
 53:   */
 54:   MatDestroy(&A);
 55:   MatDestroy(&Ar);
 56:   MatDestroy(&C);
 57:   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*/