xref: /petsc/src/mat/tests/ex4.c (revision 8fb5bd83c3955fefcf33a54e3bb66920a9fa884b)
1 
2 static char help[] = "Creates a matrix, inserts some values, and tests MatCreateSubMatrices() and MatZeroEntries().\n\n";
3 
4 #include <petscmat.h>
5 
6 int main(int argc,char **argv)
7 {
8   Mat            mat,submat,submat1,*submatrices;
9   PetscInt       m = 10,n = 10,i = 4,tmp,rstart,rend;
10   IS             irow,icol;
11   PetscScalar    value = 1.0;
12   PetscViewer    sviewer;
13   PetscBool      allA = PETSC_FALSE;
14 
15   PetscCall(PetscInitialize(&argc,&argv,(char*)0,help));
16   PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON));
17   PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_COMMON));
18 
19   PetscCall(MatCreate(PETSC_COMM_WORLD,&mat));
20   PetscCall(MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n));
21   PetscCall(MatSetFromOptions(mat));
22   PetscCall(MatSetUp(mat));
23   PetscCall(MatGetOwnershipRange(mat,&rstart,&rend));
24   for (i=rstart; i<rend; i++) {
25     value = (PetscReal)i+1; tmp = i % 5;
26     PetscCall(MatSetValues(mat,1,&tmp,1,&i,&value,INSERT_VALUES));
27   }
28   PetscCall(MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY));
29   PetscCall(MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY));
30   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original matrix\n"));
31   PetscCall(MatView(mat,PETSC_VIEWER_STDOUT_WORLD));
32 
33   /* Test MatCreateSubMatrix_XXX_All(), i.e., submatrix = A */
34   PetscCall(PetscOptionsGetBool(NULL,NULL,"-test_all",&allA,NULL));
35   if (allA) {
36     PetscCall(ISCreateStride(PETSC_COMM_SELF,m,0,1,&irow));
37     PetscCall(ISCreateStride(PETSC_COMM_SELF,n,0,1,&icol));
38     PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices));
39     PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices));
40     submat = *submatrices;
41 
42     /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */
43     PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with all\n"));
44     PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"--------------------\n"));
45     PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
46     PetscCall(MatView(submat,sviewer));
47     PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
48     PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
49 
50     PetscCall(ISDestroy(&irow));
51     PetscCall(ISDestroy(&icol));
52 
53     /* test getting a reference on a submat */
54     PetscCall(PetscObjectReference((PetscObject)submat));
55     PetscCall(MatDestroySubMatrices(1,&submatrices));
56     PetscCall(MatDestroy(&submat));
57   }
58 
59   /* Form submatrix with rows 2-4 and columns 4-8 */
60   PetscCall(ISCreateStride(PETSC_COMM_SELF,3,2,1,&irow));
61   PetscCall(ISCreateStride(PETSC_COMM_SELF,5,4,1,&icol));
62   PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices));
63   submat = *submatrices;
64 
65   /* Test reuse submatrices */
66   PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices));
67 
68   /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */
69   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices\n"));
70   PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
71   PetscCall(MatView(submat,sviewer));
72   PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
73   PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
74   PetscCall(PetscObjectReference((PetscObject)submat));
75   PetscCall(MatDestroySubMatrices(1,&submatrices));
76   PetscCall(MatDestroy(&submat));
77 
78   /* Form submatrix with rows 2-4 and all columns */
79   PetscCall(ISDestroy(&icol));
80   PetscCall(ISCreateStride(PETSC_COMM_SELF,10,0,1,&icol));
81   PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices));
82   PetscCall(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices));
83   submat = *submatrices;
84 
85   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with allcolumns\n"));
86   PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
87   PetscCall(MatView(submat,sviewer));
88   PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer));
89   PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
90 
91   /* Test MatDuplicate */
92   PetscCall(MatDuplicate(submat,MAT_COPY_VALUES,&submat1));
93   PetscCall(MatDestroy(&submat1));
94 
95   /* Zero the original matrix */
96   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original zeroed matrix\n"));
97   PetscCall(MatZeroEntries(mat));
98   PetscCall(MatView(mat,PETSC_VIEWER_STDOUT_WORLD));
99 
100   PetscCall(ISDestroy(&irow));
101   PetscCall(ISDestroy(&icol));
102   PetscCall(PetscObjectReference((PetscObject)submat));
103   PetscCall(MatDestroySubMatrices(1,&submatrices));
104   PetscCall(MatDestroy(&submat));
105   PetscCall(MatDestroy(&mat));
106   PetscCall(PetscFinalize());
107   return 0;
108 }
109 
110 /*TEST
111 
112    test:
113       args: -mat_type aij
114 
115    test:
116       suffix: 2
117       args: -mat_type dense
118 
119    test:
120       suffix: 3
121       nsize: 3
122       args: -mat_type aij
123 
124    test:
125       suffix: 4
126       nsize: 3
127       args: -mat_type dense
128 
129    test:
130       suffix: 5
131       nsize: 3
132       args: -mat_type aij -test_all
133 
134 TEST*/
135