1 static char help[] = "Test MatTransposeColoring for SeqAIJ matrices. Used for '-matmattransmult_color' on MatMatTransposeMult \n\n"; 2 3 #include <petscmat.h> 4 #include <petsc/private/matimpl.h> /* Need struct _p_MatTransposeColoring for this test. */ 5 6 int main(int argc, char **argv) { 7 Mat A, R, C, C_dense, C_sparse, Rt_dense, P, PtAP; 8 PetscInt row, col, m, n; 9 MatScalar one = 1.0, val; 10 MatColoring mc; 11 MatTransposeColoring matcoloring = 0; 12 ISColoring iscoloring; 13 PetscBool equal; 14 PetscMPIInt size; 15 16 PetscFunctionBeginUser; 17 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 18 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 19 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 20 21 /* Create seqaij A */ 22 PetscCall(MatCreate(PETSC_COMM_SELF, &A)); 23 PetscCall(MatSetSizes(A, 4, 4, 4, 4)); 24 PetscCall(MatSetType(A, MATSEQAIJ)); 25 PetscCall(MatSetFromOptions(A)); 26 PetscCall(MatSetUp(A)); 27 row = 0; 28 col = 0; 29 val = 1.0; 30 PetscCall(MatSetValues(A, 1, &row, 1, &col, &val, ADD_VALUES)); 31 row = 1; 32 col = 3; 33 val = 2.0; 34 PetscCall(MatSetValues(A, 1, &row, 1, &col, &val, ADD_VALUES)); 35 row = 2; 36 col = 2; 37 val = 3.0; 38 PetscCall(MatSetValues(A, 1, &row, 1, &col, &val, ADD_VALUES)); 39 row = 3; 40 col = 0; 41 val = 4.0; 42 PetscCall(MatSetValues(A, 1, &row, 1, &col, &val, ADD_VALUES)); 43 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 44 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 45 PetscCall(MatSetOptionsPrefix(A, "A_")); 46 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 47 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 48 49 /* Create seqaij R */ 50 PetscCall(MatCreate(PETSC_COMM_SELF, &R)); 51 PetscCall(MatSetSizes(R, 2, 4, 2, 4)); 52 PetscCall(MatSetType(R, MATSEQAIJ)); 53 PetscCall(MatSetFromOptions(R)); 54 PetscCall(MatSetUp(R)); 55 row = 0; 56 col = 0; 57 PetscCall(MatSetValues(R, 1, &row, 1, &col, &one, ADD_VALUES)); 58 row = 0; 59 col = 1; 60 PetscCall(MatSetValues(R, 1, &row, 1, &col, &one, ADD_VALUES)); 61 62 row = 1; 63 col = 1; 64 PetscCall(MatSetValues(R, 1, &row, 1, &col, &one, ADD_VALUES)); 65 row = 1; 66 col = 2; 67 PetscCall(MatSetValues(R, 1, &row, 1, &col, &one, ADD_VALUES)); 68 row = 1; 69 col = 3; 70 PetscCall(MatSetValues(R, 1, &row, 1, &col, &one, ADD_VALUES)); 71 PetscCall(MatAssemblyBegin(R, MAT_FINAL_ASSEMBLY)); 72 PetscCall(MatAssemblyEnd(R, MAT_FINAL_ASSEMBLY)); 73 PetscCall(MatSetOptionsPrefix(R, "R_")); 74 PetscCall(MatView(R, PETSC_VIEWER_STDOUT_WORLD)); 75 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 76 77 /* C = A*R^T */ 78 PetscCall(MatMatTransposeMult(A, R, MAT_INITIAL_MATRIX, 2.0, &C)); 79 PetscCall(MatSetOptionsPrefix(C, "ARt_")); 80 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD)); 81 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 82 83 /* Create MatTransposeColoring from symbolic C=A*R^T */ 84 PetscCall(MatColoringCreate(C, &mc)); 85 PetscCall(MatColoringSetDistance(mc, 2)); 86 /* PetscCall(MatColoringSetType(mc,MATCOLORINGSL)); */ 87 PetscCall(MatColoringSetFromOptions(mc)); 88 PetscCall(MatColoringApply(mc, &iscoloring)); 89 PetscCall(MatColoringDestroy(&mc)); 90 PetscCall(MatTransposeColoringCreate(C, iscoloring, &matcoloring)); 91 PetscCall(ISColoringDestroy(&iscoloring)); 92 93 /* Create Rt_dense */ 94 PetscCall(MatCreate(PETSC_COMM_WORLD, &Rt_dense)); 95 PetscCall(MatSetSizes(Rt_dense, 4, matcoloring->ncolors, PETSC_DECIDE, PETSC_DECIDE)); 96 PetscCall(MatSetType(Rt_dense, MATDENSE)); 97 PetscCall(MatSeqDenseSetPreallocation(Rt_dense, NULL)); 98 PetscCall(MatAssemblyBegin(Rt_dense, MAT_FINAL_ASSEMBLY)); 99 PetscCall(MatAssemblyEnd(Rt_dense, MAT_FINAL_ASSEMBLY)); 100 PetscCall(MatGetLocalSize(Rt_dense, &m, &n)); 101 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Rt_dense: %" PetscInt_FMT ",%" PetscInt_FMT "\n", m, n)); 102 103 /* Get Rt_dense by Apply MatTransposeColoring to R */ 104 PetscCall(MatTransColoringApplySpToDen(matcoloring, R, Rt_dense)); 105 106 /* C_dense = A*Rt_dense */ 107 PetscCall(MatMatMult(A, Rt_dense, MAT_INITIAL_MATRIX, 2.0, &C_dense)); 108 PetscCall(MatSetOptionsPrefix(C_dense, "ARt_dense_")); 109 /*PetscCall(MatView(C_dense,PETSC_VIEWER_STDOUT_WORLD)); */ 110 /*PetscCall(PetscPrintf(PETSC_COMM_SELF,"\n")); */ 111 112 /* Recover C from C_dense */ 113 PetscCall(MatDuplicate(C, MAT_DO_NOT_COPY_VALUES, &C_sparse)); 114 PetscCall(MatTransColoringApplyDenToSp(matcoloring, C_dense, C_sparse)); 115 PetscCall(MatSetOptionsPrefix(C_sparse, "ARt_color_")); 116 PetscCall(MatView(C_sparse, PETSC_VIEWER_STDOUT_WORLD)); 117 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 118 119 PetscCall(MatDestroy(&C_dense)); 120 PetscCall(MatDestroy(&C_sparse)); 121 PetscCall(MatDestroy(&Rt_dense)); 122 PetscCall(MatTransposeColoringDestroy(&matcoloring)); 123 PetscCall(MatDestroy(&C)); 124 125 /* Test PtAP = P^T*A*P, P = R^T */ 126 PetscCall(MatTranspose(R, MAT_INITIAL_MATRIX, &P)); 127 PetscCall(MatPtAP(A, P, MAT_INITIAL_MATRIX, 2.0, &PtAP)); 128 PetscCall(MatSetOptionsPrefix(PtAP, "PtAP_")); 129 PetscCall(MatView(PtAP, PETSC_VIEWER_STDOUT_WORLD)); 130 PetscCall(MatDestroy(&P)); 131 132 /* Test C = RARt */ 133 PetscCall(MatRARt(A, R, MAT_INITIAL_MATRIX, 2.0, &C)); 134 PetscCall(MatRARt(A, R, MAT_REUSE_MATRIX, 2.0, &C)); 135 PetscCall(MatEqual(C, PtAP, &equal)); 136 PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error: PtAP != RARt"); 137 138 /* Free spaces */ 139 PetscCall(MatDestroy(&C)); 140 PetscCall(MatDestroy(&A)); 141 PetscCall(MatDestroy(&R)); 142 PetscCall(MatDestroy(&PtAP)); 143 PetscCall(PetscFinalize()); 144 return 0; 145 } 146 147 /*TEST 148 149 test: 150 output_file: output/ex161.out 151 test: 152 suffix: 2 153 args: -matmattransmult_via color 154 output_file: output/ex161.out 155 156 test: 157 suffix: 3 158 args: -matmattransmult_via color -matden2sp_brows 3 159 output_file: output/ex161.out 160 161 test: 162 suffix: 4 163 args: -matmattransmult_via color -matrart_via r*art 164 output_file: output/ex161.out 165 166 test: 167 suffix: 5 168 args: -matmattransmult_via color -matrart_via coloring_rart 169 output_file: output/ex161.out 170 171 TEST*/ 172