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