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