Actual source code: ex248.c

  1: static char help[] = "Tests MatSeqAIJKron.\n\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Mat                A, B, C, K, Ad, Bd;
  8:   const PetscScalar *Bv;
  9:   PetscInt           n = 10, m = 20, p = 7, q = 17;
 10:   PetscBool          flg;

 12:   PetscFunctionBeginUser;
 13:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 14:   PetscCall(MatCreateDense(PETSC_COMM_SELF, m, n, m, n, NULL, &Ad));
 15:   PetscCall(MatCreateDense(PETSC_COMM_SELF, p, q, p, q, NULL, &Bd));
 16:   PetscCall(MatSetRandom(Ad, NULL));
 17:   PetscCall(MatSetRandom(Bd, NULL));
 18:   PetscCall(MatFilter(Ad, 0.2, PETSC_FALSE, PETSC_FALSE));
 19:   PetscCall(MatFilter(Bd, 0.2, PETSC_FALSE, PETSC_FALSE));
 20:   PetscCall(MatConvert(Ad, MATAIJ, MAT_INITIAL_MATRIX, &A));
 21:   PetscCall(MatConvert(Bd, MATAIJ, MAT_INITIAL_MATRIX, &B));
 22:   PetscCall(MatSeqAIJKron(A, B, MAT_INITIAL_MATRIX, &C));
 23:   PetscCall(MatViewFromOptions(A, NULL, "-A_view"));
 24:   PetscCall(MatViewFromOptions(B, NULL, "-B_view"));
 25:   PetscCall(MatViewFromOptions(C, NULL, "-C_view"));
 26:   PetscCall(MatDenseGetArrayRead(Bd, &Bv));
 27:   PetscCall(MatCreateKAIJ(A, p, q, NULL, Bv, &K));
 28:   PetscCall(MatDenseRestoreArrayRead(Bd, &Bv));
 29:   PetscCall(MatMultEqual(C, K, 10, &flg));
 30:   PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "K*x != C*x");
 31:   PetscCall(MatScale(A, 1.3));
 32:   PetscCall(MatScale(B, 0.3));
 33:   PetscCall(MatScale(Bd, 0.3));
 34:   PetscCall(MatSeqAIJKron(A, B, MAT_REUSE_MATRIX, &C));
 35:   PetscCall(MatDenseGetArrayRead(Bd, &Bv));
 36:   PetscCall(MatKAIJSetT(K, p, q, Bv));
 37:   PetscCall(MatDenseRestoreArrayRead(Bd, &Bv));
 38:   PetscCall(MatMultEqual(C, K, 10, &flg));
 39:   PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "K*x != C*x");
 40:   PetscCall(MatDestroy(&K));
 41:   PetscCall(MatDestroy(&C));
 42:   PetscCall(MatDestroy(&B));
 43:   PetscCall(MatDestroy(&A));
 44:   PetscCall(MatDestroy(&Bd));
 45:   PetscCall(MatDestroy(&Ad));
 46:   PetscCall(PetscFinalize());
 47:   return 0;
 48: }

 50: /*TEST

 52:     test:
 53:       suffix: 1
 54:       nsize: 1

 56: TEST*/