1 static char help[] = "Test MatProductReplaceMats() \n\ 2 Modified from the code contributed by Pierre Jolivet \n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) { 7 PetscInt n = 2, convert; 8 Mat A, B, Bdense, Conjugate; 9 PetscBool conjugate = PETSC_FALSE, equal, flg; 10 11 PetscFunctionBeginUser; 12 PetscCall(PetscInitialize(&argc, &args, NULL, help)); 13 14 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 15 PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n)); 16 PetscCall(MatSetType(A, MATDENSE)); 17 PetscCall(MatSetFromOptions(A)); 18 PetscCall(MatSeqDenseSetPreallocation(A, NULL)); 19 PetscCall(MatMPIDenseSetPreallocation(A, NULL)); 20 PetscCall(MatSetRandom(A, NULL)); 21 PetscCall(MatViewFromOptions(A, NULL, "-A_view")); 22 PetscCall(PetscOptionsGetBool(NULL, NULL, "-conjugate", &conjugate, NULL)); 23 24 for (convert = 0; convert < 2; convert++) { 25 /* convert dense matrix A to aij format */ 26 if (convert) PetscCall(MatConvert(A, MATAIJ, MAT_INPLACE_MATRIX, &A)); 27 28 /* compute B = A^T * A or B = A^H * A */ 29 PetscCall(MatProductCreate(A, A, NULL, &B)); 30 31 flg = PETSC_FALSE; 32 PetscCall(PetscOptionsGetBool(NULL, NULL, "-atb", &flg, NULL)); 33 if (flg) { 34 PetscCall(MatProductSetType(B, MATPRODUCT_AtB)); 35 } else { 36 PetscCall(PetscOptionsGetBool(NULL, NULL, "-ptap", &flg, NULL)); 37 if (flg) { 38 PetscCall(MatProductSetType(B, MATPRODUCT_PtAP)); 39 } else { 40 PetscCall(PetscOptionsGetBool(NULL, NULL, "-abt", &flg, NULL)); 41 if (flg) { 42 PetscCall(MatProductSetType(B, MATPRODUCT_ABt)); 43 } else { 44 PetscCall(MatProductSetType(B, MATPRODUCT_AB)); 45 } 46 } 47 } 48 PetscCall(MatProductSetFromOptions(B)); 49 PetscCall(MatProductSymbolic(B)); 50 51 PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &Conjugate)); 52 if (conjugate) PetscCall(MatConjugate(Conjugate)); 53 54 /* replace input A by Conjugate */ 55 PetscCall(MatProductReplaceMats(Conjugate, NULL, NULL, B)); 56 57 PetscCall(MatProductNumeric(B)); 58 PetscCall(MatViewFromOptions(B, NULL, "-product_view")); 59 60 PetscCall(MatDestroy(&Conjugate)); 61 if (!convert) { 62 Bdense = B; 63 B = NULL; 64 } 65 } 66 67 /* Compare Bdense and B */ 68 PetscCall(MatMultEqual(Bdense, B, 10, &equal)); 69 PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Bdense != B"); 70 71 PetscCall(MatDestroy(&Bdense)); 72 PetscCall(MatDestroy(&B)); 73 PetscCall(MatDestroy(&A)); 74 PetscCall(PetscFinalize()); 75 return 0; 76 } 77 78 /*TEST 79 80 test: 81 suffix: 1 82 args: -conjugate false -atb 83 output_file: output/ex258_1.out 84 85 test: 86 suffix: 2 87 args: -conjugate true -atb 88 output_file: output/ex258_1.out 89 90 test: 91 suffix: 3 92 args: -conjugate false 93 output_file: output/ex258_1.out 94 95 test: 96 suffix: 4 97 args: -ptap 98 output_file: output/ex258_1.out 99 100 test: 101 suffix: 5 102 args: -abt 103 output_file: output/ex258_1.out 104 105 test: 106 suffix: 6 107 nsize: 2 108 args: -conjugate false -atb 109 output_file: output/ex258_1.out 110 111 test: 112 suffix: 7 113 nsize: 2 114 args: -conjugate true -atb 115 output_file: output/ex258_1.out 116 117 test: 118 suffix: 8 119 nsize: 2 120 args: -conjugate false 121 output_file: output/ex258_1.out 122 123 test: 124 suffix: 9 125 nsize: 2 126 args: -ptap 127 output_file: output/ex258_1.out 128 129 test: 130 suffix: 10 131 nsize: 2 132 args: -abt 133 output_file: output/ex258_1.out 134 135 test: 136 suffix: 11 137 nsize: 2 138 args: -conjugate true -atb -mat_product_algorithm backend 139 TODO: bug: MatProductReplaceMats() does not change the product for this test 140 141 TEST*/ 142