1 2 static char help[] = "Tests MatOption MAT_FORCE_DIAGONAL_ENTRIES.\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc,char **args) 7 { 8 Mat A,B; 9 Vec diag; 10 PetscInt i,n = 10,col[3],test; 11 PetscScalar v[3]; 12 13 CHKERRQ(PetscInitialize(&argc,&args,(char*)0,help)); 14 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 15 16 /* Create A which has empty 0-th row and column */ 17 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A)); 18 CHKERRQ(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n)); 19 CHKERRQ(MatSetType(A,MATAIJ)); 20 CHKERRQ(MatSetFromOptions(A)); 21 CHKERRQ(MatSetUp(A)); 22 23 v[0] = -1.; v[1] = 2.; v[2] = -1.; 24 for (i=2; i<n-1; i++) { 25 col[0] = i-1; col[1] = i; col[2] = i+1; 26 CHKERRQ(MatSetValues(A,1,&i,3,col,v,INSERT_VALUES)); 27 } 28 i = 1; col[0] = 1; col[1] = 2; 29 CHKERRQ(MatSetValues(A,1,&i,2,col,v+1,INSERT_VALUES)); 30 i = n - 1; col[0] = n - 2; col[1] = n - 1; 31 CHKERRQ(MatSetValues(A,1,&i,2,col,v,INSERT_VALUES)); 32 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 33 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 34 35 for (test = 0; test < 2; test++) { 36 CHKERRQ(MatProductCreate(A,A,NULL,&B)); 37 38 if (test == 0) { 39 /* Compute B = A*A; B misses 0-th diagonal */ 40 CHKERRQ(MatProductSetType(B,MATPRODUCT_AB)); 41 CHKERRQ(MatSetOptionsPrefix(B,"AB_")); 42 } else { 43 /* Compute B = A^t*A; B misses 0-th diagonal */ 44 CHKERRQ(MatProductSetType(B,MATPRODUCT_AtB)); 45 CHKERRQ(MatSetOptionsPrefix(B,"AtB_")); 46 } 47 48 /* Force allocate missing diagonal entries of B */ 49 CHKERRQ(MatSetOption(B,MAT_FORCE_DIAGONAL_ENTRIES,PETSC_TRUE)); 50 CHKERRQ(MatProductSetFromOptions(B)); 51 52 CHKERRQ(MatProductSymbolic(B)); 53 CHKERRQ(MatProductNumeric(B)); 54 55 CHKERRQ(MatSetOption(B,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE)); 56 57 /* Insert entries to diagonal of B */ 58 CHKERRQ(MatCreateVecs(B,NULL,&diag)); 59 CHKERRQ(MatGetDiagonal(B,diag)); 60 CHKERRQ(VecSetValue(diag,0,100.0,INSERT_VALUES)); 61 CHKERRQ(VecAssemblyBegin(diag)); 62 CHKERRQ(VecAssemblyEnd(diag)); 63 64 CHKERRQ(MatDiagonalSet(B,diag,INSERT_VALUES)); 65 if (test == 1) { 66 CHKERRQ(MatView(B,PETSC_VIEWER_STDOUT_WORLD)); 67 } 68 CHKERRQ(MatDestroy(&B)); 69 CHKERRQ(VecDestroy(&diag)); 70 } 71 72 CHKERRQ(MatDestroy(&A)); 73 CHKERRQ(PetscFinalize()); 74 return 0; 75 } 76 77 /*TEST 78 79 test: 80 output_file: output/ex81_1.out 81 82 test: 83 suffix: 2 84 args: -AtB_mat_product_algorithm at*b 85 output_file: output/ex81_1.out 86 87 test: 88 suffix: 3 89 args: -AtB_mat_product_algorithm outerproduct 90 output_file: output/ex81_1.out 91 92 test: 93 suffix: 4 94 nsize: 3 95 args: -AtB_mat_product_algorithm nonscalable 96 output_file: output/ex81_3.out 97 98 test: 99 suffix: 5 100 nsize: 3 101 args: -AtB_mat_product_algorithm scalable 102 output_file: output/ex81_3.out 103 104 test: 105 suffix: 6 106 nsize: 3 107 args: -AtB_mat_product_algorithm at*b 108 output_file: output/ex81_3.out 109 110 TEST*/ 111