Actual source code: ex16.c

  2: static char help[] = "Tests MatGetArray() and MatView_SeqDense_Binary(), MatView_MPIDense_Binary().\n\n";

 4:  #include petscmat.h
 5:  #include petscviewer.h

  9: int main(int argc,char **args)
 10: {
 11:   Mat              A;
 12:   PetscInt         i,j,m = 3,n = 2,rstart,rend;
 13:   PetscErrorCode   ierr;
 14:   PetscScalar      v,*array;
 15:   PetscViewer      view;

 17:   PetscInitialize(&argc,&args,(char *)0,help);

 19:   /*
 20:       Create a parallel dense matrix shared by all processors 
 21:   */
 22:   MatCreateMPIDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,PETSC_NULL,&A);
 23: 

 25:   /*
 26:      Set values into the matrix 
 27:   */
 28:   for (i=0; i<m; i++) {
 29:     for (j=0; j<n; j++) {
 30:       v = 9.0/(i+j+1); MatSetValues(A,1,&i,1,&j,&v,INSERT_VALUES);
 31:     }
 32:   }
 33:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 34:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 36:   /*
 37:        Print the matrix to the screen 
 38:   */
 39:   MatView(A,PETSC_VIEWER_STDOUT_WORLD);


 42:   /*
 43:       Print the local portion of the matrix to the screen
 44:   */
 45:   MatGetArray(A,&array);
 46:   MatGetOwnershipRange(A,&rstart,&rend);
 47:   for (i=rstart; i<rend; i++) {
 48:     for (j=0; j<n; j++) {
 49:       PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%6.4e ",PetscRealPart(array[j*(rend-rstart)+i-rstart]));
 50:     }
 51:     PetscSynchronizedPrintf(PETSC_COMM_WORLD,"\n");
 52:   }
 53:   PetscSynchronizedFlush(PETSC_COMM_WORLD);
 54:   MatRestoreArray(A,&array);

 56:   /*
 57:       Store the binary matrix to a file
 58:   */
 59:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);
 60:   PetscViewerSetFormat(view, PETSC_VIEWER_NATIVE);
 61:   MatView(A,view);
 62:   PetscViewerDestroy(view);
 63:   MatDestroy(A);

 65:   /* 
 66:      Now reload the matrix and view it
 67:   */
 68:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);
 69:   MatLoad(view,MATMPIDENSE,&A);
 70:   PetscViewerDestroy(view);
 71:   MatView(A,PETSC_VIEWER_STDOUT_WORLD);
 72:   MatDestroy(A);

 74:   PetscFinalize();
 75:   return 0;
 76: }