Actual source code: ex248.c


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

  4: #include <petscmat.h>

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

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

 51: /*TEST

 53:     test:
 54:       suffix: 1
 55:       nsize: 1
 56:       output_file: output/ex101.out

 58: TEST*/