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 { 8 PetscInt n = 2,convert; 9 Mat A,B,Bdense,Conjugate; 10 PetscBool conjugate = PETSC_FALSE,equal,flg; 11 12 PetscFunctionBeginUser; 13 PetscCall(PetscInitialize(&argc,&args,NULL,help)); 14 15 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 16 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n)); 17 PetscCall(MatSetType(A,MATDENSE)); 18 PetscCall(MatSetFromOptions(A)); 19 PetscCall(MatSeqDenseSetPreallocation(A,NULL)); 20 PetscCall(MatMPIDenseSetPreallocation(A,NULL)); 21 PetscCall(MatSetRandom(A,NULL)); 22 PetscCall(MatViewFromOptions(A,NULL,"-A_view")); 23 PetscCall(PetscOptionsGetBool(NULL,NULL,"-conjugate",&conjugate,NULL)); 24 25 for (convert = 0; convert<2; convert++) { 26 /* convert dense matrix A to aij format */ 27 if (convert) PetscCall(MatConvert(A,MATAIJ,MAT_INPLACE_MATRIX,&A)); 28 29 /* compute B = A^T * A or B = A^H * A */ 30 PetscCall(MatProductCreate(A,A,NULL,&B)); 31 32 flg = PETSC_FALSE; 33 PetscCall(PetscOptionsGetBool(NULL,NULL,"-atb",&flg,NULL)); 34 if (flg) { 35 PetscCall(MatProductSetType(B,MATPRODUCT_AtB)); 36 } else { 37 PetscCall(PetscOptionsGetBool(NULL,NULL,"-ptap",&flg,NULL)); 38 if (flg) { 39 PetscCall(MatProductSetType(B,MATPRODUCT_PtAP)); 40 } else { 41 PetscCall(PetscOptionsGetBool(NULL,NULL,"-abt",&flg,NULL)); 42 if (flg) { 43 PetscCall(MatProductSetType(B,MATPRODUCT_ABt)); 44 } else { 45 PetscCall(MatProductSetType(B,MATPRODUCT_AB)); 46 } 47 } 48 } 49 PetscCall(MatProductSetFromOptions(B)); 50 PetscCall(MatProductSymbolic(B)); 51 52 PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &Conjugate)); 53 if (conjugate) PetscCall(MatConjugate(Conjugate)); 54 55 /* replace input A by Conjugate */ 56 PetscCall(MatProductReplaceMats(Conjugate,NULL,NULL,B)); 57 58 PetscCall(MatProductNumeric(B)); 59 PetscCall(MatViewFromOptions(B,NULL,"-product_view")); 60 61 PetscCall(MatDestroy(&Conjugate)); 62 if (!convert) { 63 Bdense = B; 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