Actual source code: ex68.c
1: /*$Id: ex68.c,v 1.17 2001/08/07 03:03:07 balay Exp $*/
3: static char help[] = "Tests MatReorderForNonzeroDiagonal().nn";
5: #include petscmat.h
7: int main(int argc,char **argv)
8: {
9: Mat mat,B;
10: int ierr,i,j;
11: PetscScalar v;
12: IS isrow,iscol;
14: PetscInitialize(&argc,&argv,(char*)0,help);
17: /* ------- Assemble matrix, --------- */
19: MatCreateSeqAIJ(PETSC_COMM_WORLD,4,4,0,0,&mat);
21: /* set anti-diagonal of matrix */
22: v = 1.0;
23: i = 0; j = 3;
24: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
25: v = 2.0;
26: i = 1; j = 2;
27: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
28: v = 3.0;
29: i = 2; j = 1;
30: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
31: v = 4.0;
32: i = 3; j = 0;
33: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
35: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
36: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
38: printf("Original matrixn");
39: PetscViewerSetFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_DENSE);
40: MatView(mat,PETSC_VIEWER_STDOUT_SELF);
42: MatGetOrdering(mat,MATORDERING_NATURAL,&isrow,&iscol);
44: MatPermute(mat,isrow,iscol,&B);
45: printf("Original matrix permuted by identityn");
46: MatView(B,PETSC_VIEWER_STDOUT_SELF);
47: MatDestroy(B);
49: MatReorderForNonzeroDiagonal(mat,1.e-8,isrow,iscol);
50: MatPermute(mat,isrow,iscol,&B);
51: printf("Original matrix permuted by identity + NonzeroDiagonal()n");
52: MatView(B,PETSC_VIEWER_STDOUT_SELF);
53: printf("Row permutationn");
54: ISView(isrow,PETSC_VIEWER_STDOUT_SELF);
55: printf("Column permutationn");
56: ISView(iscol,PETSC_VIEWER_STDOUT_SELF);
57: MatDestroy(B);
59: ISDestroy(isrow);
60: ISDestroy(iscol);
62: MatGetOrdering(mat,MATORDERING_ND,&isrow,&iscol);
63: MatPermute(mat,isrow,iscol,&B);
64: printf("Original matrix permuted by NDn");
65: MatView(B,PETSC_VIEWER_STDOUT_SELF);
66: MatDestroy(B);
67: printf("ND row permutationn");
68: ISView(isrow,PETSC_VIEWER_STDOUT_SELF);
69: printf("ND column permutationn");
70: ISView(iscol,PETSC_VIEWER_STDOUT_SELF);
72: MatReorderForNonzeroDiagonal(mat,1.e-8,isrow,iscol);
73: MatPermute(mat,isrow,iscol,&B);
74: printf("Original matrix permuted by ND + NonzeroDiagonal()n");
75: MatView(B,PETSC_VIEWER_STDOUT_SELF);
76: MatDestroy(B);
77: printf("ND + NonzeroDiagonal() row permutationn");
78: ISView(isrow,PETSC_VIEWER_STDOUT_SELF);
79: printf("ND + NonzeroDiagonal() column permutationn");
80: ISView(iscol,PETSC_VIEWER_STDOUT_SELF);
82: ISDestroy(isrow);
83: ISDestroy(iscol);
85: MatGetOrdering(mat,MATORDERING_RCM,&isrow,&iscol);
86: MatPermute(mat,isrow,iscol,&B);
87: printf("Original matrix permuted by RCMn");
88: MatView(B,PETSC_VIEWER_STDOUT_SELF);
89: MatDestroy(B);
90: printf("RCM row permutationn");
91: ISView(isrow,PETSC_VIEWER_STDOUT_SELF);
92: printf("RCM column permutationn");
93: ISView(iscol,PETSC_VIEWER_STDOUT_SELF);
95: MatReorderForNonzeroDiagonal(mat,1.e-8,isrow,iscol);
96: MatPermute(mat,isrow,iscol,&B);
97: printf("Original matrix permuted by RCM + NonzeroDiagonal()n");
98: MatView(B,PETSC_VIEWER_STDOUT_SELF);
99: MatDestroy(B);
100: printf("RCM + NonzeroDiagonal() row permutationn");
101: ISView(isrow,PETSC_VIEWER_STDOUT_SELF);
102: printf("RCM + NonzeroDiagonal() column permutationn");
103: ISView(iscol,PETSC_VIEWER_STDOUT_SELF);
105: MatLUFactor(mat,isrow,iscol,PETSC_NULL);
106: printf("Factored matrix permuted by RCM + NonzeroDiagonal()n");
107: MatView(mat,PETSC_VIEWER_STDOUT_SELF);
109: /* Free data structures */
110: ISDestroy(isrow);
111: ISDestroy(iscol);
112: MatDestroy(mat);
114: PetscFinalize();
115: return 0;
116: }
117: