2: static char help[] = "Tests sequential and parallel MatMatMult() and MatAXPY(...,SUBSET_NONZERO_PATTERN) \n\
  3: Input arguments are:\n\
  4:   -f <input_file>  : file to load\n\n";
  5: /* e.g., mpiexec -n 3 ./ex113 -f <file> */

  7:  #include <petscmat.h>

  9: int main(int argc,char **args)
 10: {
 11:   Mat            A,A1,A2,Mtmp,dstMat;
 12:   PetscViewer    viewer;
 14:   PetscReal      fill=4.0;
 15:   char           file[128];
 16:   PetscBool      flg;

 18:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 19:   /*  Load the matrix A */
 20:   PetscOptionsGetString(NULL,NULL,"-f",file,128,&flg);
 21:   if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate a file name for matrix A with the -f option.");

 23:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);
 24:   MatCreate(PETSC_COMM_WORLD,&A);
 25:   MatLoad(A,viewer);
 26:   PetscViewerDestroy(&viewer);

 28:   MatDuplicate(A,MAT_COPY_VALUES,&A1);
 29:   MatDuplicate(A,MAT_COPY_VALUES,&A2);

 31:   /* dstMat = A*A1*A2 */
 32:   MatMatMult(A1,A2,MAT_INITIAL_MATRIX,fill,&Mtmp);
 33:   MatMatMult(A,Mtmp,MAT_INITIAL_MATRIX,fill,&dstMat);
 34:   MatDestroy(&Mtmp);

 36:   /* dstMat += A1*A2 */
 37:   MatMatMult(A1,A2,MAT_INITIAL_MATRIX,fill,&Mtmp);
 38:   MatAXPY(dstMat,1.0,Mtmp,SUBSET_NONZERO_PATTERN);
 39:   MatDestroy(&Mtmp);

 41:   /* dstMat += A*A1 */
 42:   MatMatMult(A,A1,MAT_INITIAL_MATRIX,fill,&Mtmp);
 43:   MatAXPY(dstMat, 1.0, Mtmp,SUBSET_NONZERO_PATTERN);
 44:   MatDestroy(&Mtmp);

 46:   /* dstMat += A */
 47:   MatAXPY(dstMat, 1.0, A,SUBSET_NONZERO_PATTERN);

 49:   MatDestroy(&A);
 50:   MatDestroy(&A1);
 51:   MatDestroy(&A2);
 52:   MatDestroy(&dstMat);
 53:   PetscFinalize();
 54:   return ierr;
 55: }