117ea310bSPierre Jolivet #include <petsc.h> 217ea310bSPierre Jolivet 317ea310bSPierre Jolivet static char help[] = "Tests for MatEliminateZeros().\n\n"; 417ea310bSPierre Jolivet 517ea310bSPierre Jolivet int main(int argc, char **args) 617ea310bSPierre Jolivet { 7*2ce66baaSPierre Jolivet Mat A, B, C, D, E; 817ea310bSPierre Jolivet PetscInt M = 40, bs = 2; 917ea310bSPierre Jolivet PetscReal threshold = 1.2; 1017ea310bSPierre Jolivet PetscBool flg; 1117ea310bSPierre Jolivet 1217ea310bSPierre Jolivet PetscFunctionBeginUser; 1317ea310bSPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help)); 1417ea310bSPierre Jolivet PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL)); 1517ea310bSPierre Jolivet PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL)); 1617ea310bSPierre Jolivet PetscCall(PetscOptionsGetReal(NULL, NULL, "-threshold", &threshold, NULL)); 1717ea310bSPierre Jolivet PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, M, NULL, &D)); 1817ea310bSPierre Jolivet PetscCall(MatSetRandom(D, NULL)); 1917ea310bSPierre Jolivet PetscCall(MatTranspose(D, MAT_INITIAL_MATRIX, &A)); 2017ea310bSPierre Jolivet PetscCall(MatAXPY(D, 1.0, A, SAME_NONZERO_PATTERN)); 2117ea310bSPierre Jolivet PetscCall(MatDestroy(&A)); 2217ea310bSPierre Jolivet PetscCall(MatSetBlockSize(D, bs)); 23*2ce66baaSPierre Jolivet PetscCall(MatDuplicate(D, MAT_COPY_VALUES, &E)); 24*2ce66baaSPierre Jolivet PetscCall(MatViewFromOptions(D, NULL, "-input_dense")); 25*2ce66baaSPierre Jolivet for (PetscInt i = 0; i < 2; ++i) { 2617ea310bSPierre Jolivet PetscCall(MatConvert(D, MATAIJ, MAT_INITIAL_MATRIX, &A)); 2717ea310bSPierre Jolivet PetscCall(MatConvert(D, MATBAIJ, MAT_INITIAL_MATRIX, &B)); 2817ea310bSPierre Jolivet PetscCall(MatConvert(D, MATSBAIJ, MAT_INITIAL_MATRIX, &C)); 29*2ce66baaSPierre Jolivet if (i == 0) { // filtering, elimination, but no compression 3017ea310bSPierre Jolivet PetscCall(MatViewFromOptions(A, NULL, "-input_aij")); 3117ea310bSPierre Jolivet PetscCall(MatViewFromOptions(B, NULL, "-input_baij")); 3217ea310bSPierre Jolivet PetscCall(MatViewFromOptions(C, NULL, "-input_sbaij")); 33*2ce66baaSPierre Jolivet PetscCall(MatFilter(D, threshold, PETSC_FALSE, PETSC_FALSE)); 34*2ce66baaSPierre Jolivet PetscCall(MatFilter(A, threshold, PETSC_FALSE, PETSC_FALSE)); 3517ea310bSPierre Jolivet PetscCall(MatEliminateZeros(A, PETSC_TRUE)); 36*2ce66baaSPierre Jolivet PetscCall(MatFilter(B, threshold, PETSC_FALSE, PETSC_FALSE)); 3717ea310bSPierre Jolivet PetscCall(MatEliminateZeros(B, PETSC_TRUE)); 38*2ce66baaSPierre Jolivet PetscCall(MatFilter(C, threshold, PETSC_FALSE, PETSC_FALSE)); 3917ea310bSPierre Jolivet PetscCall(MatEliminateZeros(C, PETSC_TRUE)); 40*2ce66baaSPierre Jolivet } else { // filtering, elimination, and compression 41*2ce66baaSPierre Jolivet PetscCall(MatFilter(D, threshold, PETSC_TRUE, PETSC_FALSE)); 42*2ce66baaSPierre Jolivet PetscCall(MatFilter(A, threshold, PETSC_TRUE, PETSC_FALSE)); 43*2ce66baaSPierre Jolivet PetscCall(MatFilter(B, threshold, PETSC_TRUE, PETSC_FALSE)); 44*2ce66baaSPierre Jolivet PetscCall(MatFilter(C, threshold, PETSC_TRUE, PETSC_FALSE)); 45*2ce66baaSPierre Jolivet } 4617ea310bSPierre Jolivet PetscCall(MatViewFromOptions(D, NULL, "-output_dense")); 4717ea310bSPierre Jolivet PetscCall(MatViewFromOptions(A, NULL, "-output_aij")); 4817ea310bSPierre Jolivet PetscCall(MatViewFromOptions(B, NULL, "-output_baij")); 4917ea310bSPierre Jolivet PetscCall(MatViewFromOptions(C, NULL, "-output_sbaij")); 5017ea310bSPierre Jolivet PetscCall(MatMultEqual(D, A, 10, &flg)); 5117ea310bSPierre Jolivet PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "A != D"); 5217ea310bSPierre Jolivet PetscCall(MatMultEqual(D, B, 10, &flg)); 5317ea310bSPierre Jolivet PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "B != D"); 5417ea310bSPierre Jolivet PetscCall(MatMultEqual(D, C, 10, &flg)); 5517ea310bSPierre Jolivet PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "C != D"); 5617ea310bSPierre Jolivet PetscCall(MatDestroy(&C)); 5717ea310bSPierre Jolivet PetscCall(MatDestroy(&B)); 5817ea310bSPierre Jolivet PetscCall(MatDestroy(&A)); 5917ea310bSPierre Jolivet PetscCall(MatDestroy(&D)); 60*2ce66baaSPierre Jolivet D = E; 61*2ce66baaSPierre Jolivet } 6217ea310bSPierre Jolivet PetscCall(PetscFinalize()); 6317ea310bSPierre Jolivet return 0; 6417ea310bSPierre Jolivet } 6517ea310bSPierre Jolivet 6617ea310bSPierre Jolivet /*TEST 6717ea310bSPierre Jolivet 6817ea310bSPierre Jolivet test: 6917ea310bSPierre Jolivet nsize: {{1 2}} 7017ea310bSPierre Jolivet output_file: output/empty.out 7117ea310bSPierre Jolivet 7217ea310bSPierre Jolivet TEST*/ 73