Actual source code: ex113.c


  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;
 13:   PetscReal   fill = 4.0;
 14:   char        file[128];
 15:   PetscBool   flg;

 17:   PetscFunctionBeginUser;
 18:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 19:   /*  Load the matrix A */
 20:   PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
 21:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate a file name for matrix A with the -f option.");

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

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

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

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

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

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

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