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