xref: /petsc/src/mat/tests/ex114.c (revision f07e67ed86ede464a9faddc80617894c11508dd6)
1 
2 static char help[] = "Tests MatGetRowMax(), MatGetRowMin(), MatGetRowMaxAbs()\n";
3 
4 #include <petscmat.h>
5 
6 #define M 5
7 #define N 6
8 
9 int main(int argc,char **args)
10 {
11   Mat            A;
12   Vec            min,max,maxabs,e;
13   PetscInt       m,n,j,col,cols[N];
14   PetscInt       imin[M],imax[M],imaxabs[M],indices[N],row,testcase=0;
15   PetscScalar    values[N];
16   PetscErrorCode ierr;
17   MatType        type;
18   PetscMPIInt    size,rank;
19   PetscBool      doTest=PETSC_TRUE;
20   PetscReal      enorm;
21 
22   ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
23   ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
24   ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
25   ierr = PetscOptionsGetInt(NULL,NULL,"-testcase",&testcase,NULL);CHKERRQ(ierr);
26 
27   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
28   if (testcase == 1) { /* proc[0] holds entire A and other processes have no entry */
29     if (!rank) {
30       ierr = MatSetSizes(A,M,N,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
31     } else {
32       ierr = MatSetSizes(A,0,0,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
33     }
34     testcase = 0;
35   } else {
36     ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr);
37   }
38   ierr = MatSetFromOptions(A);CHKERRQ(ierr);
39   ierr = MatSetUp(A);CHKERRQ(ierr);
40 
41   if (!rank) { /* proc[0] sets matrix A */
42     for (j=0; j<N; j++) indices[j] = j;
43     switch (testcase) {
44     case 1: /* see testcast 0 */
45       break;
46     case 2:
47       row = 0;
48       values[0]  = -2.0; values[1] = -2.0; values[2] = -2.0; values[3] = -4.0; values[4] = 1.0; values[5] = 1.0;
49       ierr = MatSetValues(A,1,&row,N,indices,values,INSERT_VALUES);CHKERRQ(ierr);
50       row = 2;
51       indices[0] = 0;    indices[1] = 3;    indices[2] = 5;
52       values[0]  = -2.0; values[1]  = -2.0; values[2]  = -2.0;
53       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
54       row = 3;
55       indices[0] = 0;    indices[1] = 1;    indices[2] = 4;
56       values[0]  = -2.0; values[1]  = -2.0; values[2]  = -2.0;
57       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
58       row = 4;
59       indices[0] = 0;    indices[1] = 1;    indices[2] = 2;
60       values[0]  = -2.0; values[1]  = -2.0; values[2]  = -2.0;
61       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
62       break;
63     case 3:
64       row = 0;
65       values[0]  = -2.0; values[1] = -2.0; values[2] = -2.0;
66       ierr = MatSetValues(A,1,&row,3,indices+1,values,INSERT_VALUES);CHKERRQ(ierr);
67       row = 1;
68       values[0]  = -2.0; values[1] = -2.0; values[2] = -2.0;
69       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
70       row = 2;
71       values[0]  = -2.0; values[1] = -2.0; values[2]  = -2.0;
72       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
73       row = 3;
74       values[0]  = -2.0; values[1] = -2.0; values[2]  = -2.0; values[3] = -1.0;
75       ierr = MatSetValues(A,1,&row,4,indices,values,INSERT_VALUES);CHKERRQ(ierr);
76       row = 4;
77       values[0]  = -2.0; values[1] = -2.0; values[2]  = -2.0; values[3] = -1.0;
78       ierr = MatSetValues(A,1,&row,4,indices,values,INSERT_VALUES);CHKERRQ(ierr);
79       break;
80 
81     default:
82       row  = 0;
83       values[0]  = -1.0; values[1] = 0.0; values[2] = 1.0; values[3] = 3.0; values[4] = 4.0; values[5] = -5.0;
84       ierr = MatSetValues(A,1,&row,N,indices,values,INSERT_VALUES);CHKERRQ(ierr);
85       row  = 1;
86       ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr);
87       row  = 3;
88       ierr = MatSetValues(A,1,&row,1,indices+4,values+4,INSERT_VALUES);CHKERRQ(ierr);
89       row  = 4;
90       ierr = MatSetValues(A,1,&row,2,indices+4,values+4,INSERT_VALUES);CHKERRQ(ierr);
91     }
92   }
93   ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
94   ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
95   ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
96 
97   ierr = MatGetLocalSize(A, &m,&n);CHKERRQ(ierr);
98   ierr = VecCreate(PETSC_COMM_WORLD,&min);CHKERRQ(ierr);
99   ierr = VecSetSizes(min,m,PETSC_DECIDE);CHKERRQ(ierr);
100   ierr = VecSetFromOptions(min);CHKERRQ(ierr);
101   ierr = VecDuplicate(min,&max);CHKERRQ(ierr);
102   ierr = VecDuplicate(min,&maxabs);CHKERRQ(ierr);
103   ierr = VecDuplicate(min,&e);CHKERRQ(ierr);
104 
105   /* MatGetRowMax() */
106   ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMax\n");CHKERRQ(ierr);
107   ierr = MatGetRowMax(A,max,NULL);CHKERRQ(ierr);
108   ierr = MatGetRowMax(A,max,imax);CHKERRQ(ierr);
109   ierr = VecView(max,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
110   ierr = VecGetLocalSize(max,&n);CHKERRQ(ierr);
111   ierr = PetscIntView(n,imax,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
112 
113   /* MatGetRowMin() */
114   ierr = MatScale(A,-1.0);CHKERRQ(ierr);
115   ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
116   ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMin\n");CHKERRQ(ierr);
117   ierr = MatGetRowMin(A,min,NULL);CHKERRQ(ierr);
118   ierr = MatGetRowMin(A,min,imin);CHKERRQ(ierr);
119 
120   ierr = VecWAXPY(e,1.0,max,min);CHKERRQ(ierr); /* e = max + min */
121   ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr);
122   if (enorm > PETSC_MACHINE_EPSILON) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"max+min > PETSC_MACHINE_EPSILON ");
123   for (j = 0; j < n; j++) {
124     if (imin[j] != imax[j]) {
125       SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"imin %D != imax %D",imin[j],imax[j]);
126     }
127   }
128 
129   /* MatGetRowMaxAbs() */
130   ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMaxAbs\n");CHKERRQ(ierr);
131   ierr = MatGetRowMaxAbs(A,maxabs,NULL);CHKERRQ(ierr);
132   ierr = MatGetRowMaxAbs(A,maxabs,imaxabs);CHKERRQ(ierr);
133   ierr = VecView(maxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
134   ierr = PetscIntView(n,imaxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
135 
136   /* MatGetRowMinAbs() */
137   ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMinAbs\n");CHKERRQ(ierr);
138   ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
139   ierr = MatGetRowMinAbs(A,min,NULL);CHKERRQ(ierr);
140   ierr = MatGetRowMinAbs(A,min,imin);CHKERRQ(ierr);
141   ierr = VecView(min,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
142   ierr = PetscIntView(n,imin,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
143 
144   if (size == 1) {
145     /* Test MatGetRowMax, MatGetRowMin and MatGetRowMaxAbs for SeqDense and MPIBAIJ matrix */
146     ierr = MatConvert(A,MATDENSE,MAT_INPLACE_MATRIX,&A);CHKERRQ(ierr);
147     ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMax for seqdense matrix\n");CHKERRQ(ierr);
148     ierr = MatGetRowMax(A,max,imax);CHKERRQ(ierr);
149     ierr = VecView(max,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
150     ierr = PetscIntView(n,imax,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
151 
152     ierr = MatScale(A,-1.0);CHKERRQ(ierr);
153     ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMin for seqdense matrix\n");CHKERRQ(ierr);
154     ierr = MatGetRowMin(A,min,imin);CHKERRQ(ierr);
155     ierr = VecView(min,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
156     ierr = PetscIntView(n,imin,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
157 
158     ierr = VecWAXPY(e,1.0,max,min);CHKERRQ(ierr); /* e = max + min */
159     ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr);
160     if (enorm > PETSC_MACHINE_EPSILON) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"max+min > PETSC_MACHINE_EPSILON ");
161     for (j = 0; j < n; j++) {
162       if (imin[j] != imax[j]) {
163         SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"imin %D != imax %D for seqdense matrix",imin[j],imax[j]);
164       }
165     }
166 
167     ierr = MatGetRowMaxAbs(A,maxabs,imaxabs);CHKERRQ(ierr);
168     ierr = VecView(maxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
169     ierr = PetscIntView(n,imaxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
170 
171   } else {
172     #if 0
173     /* BAIJ */
174     ierr = PetscObjectTypeCompare((PetscObject)A,MATMPIBAIJ,&doTest);CHKERRQ(ierr);
175     if (doTest) {
176       ierr = MatGetRowMaxAbs(A,maxabs,NULL);CHKERRQ(ierr);
177       ierr = MatGetRowMaxAbs(A,maxabs,imaxabs);CHKERRQ(ierr);
178       ierr = PetscPrintf(PETSC_COMM_WORLD,"Row Maximum Absolute Values:\n");CHKERRQ(ierr);
179       ierr = VecView(maxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
180     }
181     #endif
182   }
183 
184   ierr = VecDestroy(&min);CHKERRQ(ierr);
185   ierr = VecDestroy(&max);CHKERRQ(ierr);
186   ierr = VecDestroy(&maxabs);CHKERRQ(ierr);
187   ierr = VecDestroy(&e);CHKERRQ(ierr);
188   ierr = MatDestroy(&A);CHKERRQ(ierr);
189   ierr = PetscFinalize();
190   return ierr;
191 }
192 
193 
194 
195 /*TEST
196 
197    test:
198       output_file: output/ex114.out
199 
200    test:
201       suffix: 2
202       nsize: 2
203 
204    test:
205       suffix: 3
206       nsize: 2
207       args: -mat_type baij
208 
209 TEST*/
210