Actual source code: ex216.c

  1: #include <petsc.h>

  3: static char help[] = "Tests for MatEliminateZeros().\n\n";

  5: int main(int argc, char **args)
  6: {
  7:   Mat       A, B, C, D, E;
  8:   PetscInt  M = 40, bs = 2;
  9:   PetscReal threshold = 1.2;
 10:   PetscBool flg;

 12:   PetscFunctionBeginUser;
 13:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 14:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
 15:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL));
 16:   PetscCall(PetscOptionsGetReal(NULL, NULL, "-threshold", &threshold, NULL));
 17:   PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, M, NULL, &D));
 18:   PetscCall(MatSetRandom(D, NULL));
 19:   PetscCall(MatTranspose(D, MAT_INITIAL_MATRIX, &A));
 20:   PetscCall(MatAXPY(D, 1.0, A, SAME_NONZERO_PATTERN));
 21:   PetscCall(MatDestroy(&A));
 22:   PetscCall(MatSetBlockSize(D, bs));
 23:   PetscCall(MatDuplicate(D, MAT_COPY_VALUES, &E));
 24:   PetscCall(MatViewFromOptions(D, NULL, "-input_dense"));
 25:   for (PetscInt i = 0; i < 2; ++i) {
 26:     PetscCall(MatConvert(D, MATAIJ, MAT_INITIAL_MATRIX, &A));
 27:     PetscCall(MatConvert(D, MATBAIJ, MAT_INITIAL_MATRIX, &B));
 28:     PetscCall(MatConvert(D, MATSBAIJ, MAT_INITIAL_MATRIX, &C));
 29:     if (i == 0) { // filtering, elimination, but no compression
 30:       PetscCall(MatViewFromOptions(A, NULL, "-input_aij"));
 31:       PetscCall(MatViewFromOptions(B, NULL, "-input_baij"));
 32:       PetscCall(MatViewFromOptions(C, NULL, "-input_sbaij"));
 33:       PetscCall(MatFilter(D, threshold, PETSC_FALSE, PETSC_FALSE));
 34:       PetscCall(MatFilter(A, threshold, PETSC_FALSE, PETSC_FALSE));
 35:       PetscCall(MatEliminateZeros(A, PETSC_TRUE));
 36:       PetscCall(MatFilter(B, threshold, PETSC_FALSE, PETSC_FALSE));
 37:       PetscCall(MatEliminateZeros(B, PETSC_TRUE));
 38:       PetscCall(MatFilter(C, threshold, PETSC_FALSE, PETSC_FALSE));
 39:       PetscCall(MatEliminateZeros(C, PETSC_TRUE));
 40:     } else { // filtering, elimination, and compression
 41:       PetscCall(MatFilter(D, threshold, PETSC_TRUE, PETSC_FALSE));
 42:       PetscCall(MatFilter(A, threshold, PETSC_TRUE, PETSC_FALSE));
 43:       PetscCall(MatFilter(B, threshold, PETSC_TRUE, PETSC_FALSE));
 44:       PetscCall(MatFilter(C, threshold, PETSC_TRUE, PETSC_FALSE));
 45:     }
 46:     PetscCall(MatViewFromOptions(D, NULL, "-output_dense"));
 47:     PetscCall(MatViewFromOptions(A, NULL, "-output_aij"));
 48:     PetscCall(MatViewFromOptions(B, NULL, "-output_baij"));
 49:     PetscCall(MatViewFromOptions(C, NULL, "-output_sbaij"));
 50:     PetscCall(MatMultEqual(D, A, 10, &flg));
 51:     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "A != D");
 52:     PetscCall(MatMultEqual(D, B, 10, &flg));
 53:     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "B != D");
 54:     PetscCall(MatMultEqual(D, C, 10, &flg));
 55:     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "C != D");
 56:     PetscCall(MatDestroy(&C));
 57:     PetscCall(MatDestroy(&B));
 58:     PetscCall(MatDestroy(&A));
 59:     PetscCall(MatDestroy(&D));
 60:     D = E;
 61:   }
 62:   PetscCall(PetscFinalize());
 63:   return 0;
 64: }

 66: /*TEST

 68:    test:
 69:       nsize: {{1 2}}
 70:       output_file: output/empty.out

 72: TEST*/