Actual source code: ex95.c

  1: static char help[] = "Testing MatCreateMPIAIJSumSeqAIJ().\n\n";

  3: #include <petscmat.h>

  5: int main(int argc,char **argv)
  6: {
  7:   Mat            A,B;
  8:   MatScalar      a[1],alpha;
  9:   PetscMPIInt    size,rank;
 10:   PetscInt       m,n,i,col, prid;

 12:   PetscInitialize(&argc,&argv,(char*)0,help);
 13:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 14:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 15:   prid = size;
 16:   PetscOptionsGetInt(NULL,NULL,"-prid",&prid,NULL);

 18:   m    = n = 10*size;
 19:   MatCreate(PETSC_COMM_SELF,&A);
 20:   MatSetSizes(A,PETSC_DETERMINE,PETSC_DETERMINE,m,n);
 21:   MatSetType(A,MATSEQAIJ);
 22:   MatSetUp(A);

 24:   a[0] = rank+1;
 25:   for (i=0; i<m-rank; i++) {
 26:     col  = i+rank;
 27:     MatSetValues(A,1,&i,1,&col,a,INSERT_VALUES);
 28:   }
 29:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 30:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 32:   if (rank == prid) {
 33:     PetscPrintf(PETSC_COMM_SELF,"[%d] A: \n",rank);
 34:     MatView(A,PETSC_VIEWER_STDOUT_SELF);
 35:   }

 37:   /* Test MatCreateMPIAIJSumSeqAIJ */
 38:   MatCreateMPIAIJSumSeqAIJ(PETSC_COMM_WORLD,A,PETSC_DECIDE,PETSC_DECIDE,MAT_INITIAL_MATRIX,&B);

 40:   /* Test MAT_REUSE_MATRIX */
 41:   alpha = 0.1;
 42:   for (i=0; i<3; i++) {
 43:     MatScale(A,alpha);
 44:     MatCreateMPIAIJSumSeqAIJ(PETSC_COMM_WORLD,A,PETSC_DECIDE,PETSC_DECIDE,MAT_REUSE_MATRIX,&B);
 45:   }
 46:   MatView(B, PETSC_VIEWER_STDOUT_WORLD);
 47:   MatDestroy(&B);
 48:   MatDestroy(&A);
 49:   PetscFinalize();
 50:   return 0;
 51: }

 53: /*TEST

 55:    test:
 56:       nsize: 3
 57:       filter: grep -v "MPI processes"

 59:    test:
 60:       suffix: 2
 61:       filter: grep -v "MPI processes"

 63: TEST*/