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,imin[M],imax[M],imaxabs[M],indices[N],row,testcase=0; 14 PetscScalar values[N]; 15 PetscErrorCode ierr; 16 PetscMPIInt size,rank; 17 PetscReal enorm; 18 19 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 20 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 21 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 22 ierr = PetscOptionsGetInt(NULL,NULL,"-testcase",&testcase,NULL);CHKERRQ(ierr); 23 24 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 25 if (testcase == 1) { /* proc[0] holds entire A and other processes have no entry */ 26 if (!rank) { 27 ierr = MatSetSizes(A,M,N,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr); 28 } else { 29 ierr = MatSetSizes(A,0,0,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr); 30 } 31 testcase = 0; 32 } else { 33 ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr); 34 } 35 ierr = MatSetFromOptions(A);CHKERRQ(ierr); 36 ierr = MatSetUp(A);CHKERRQ(ierr); 37 38 if (!rank) { /* proc[0] sets matrix A */ 39 for (j=0; j<N; j++) indices[j] = j; 40 switch (testcase) { 41 case 1: /* see testcast 0 */ 42 break; 43 case 2: 44 row = 0; 45 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; values[3] = -4.0; values[4] = 1.0; values[5] = 1.0; 46 ierr = MatSetValues(A,1,&row,N,indices,values,INSERT_VALUES);CHKERRQ(ierr); 47 row = 2; 48 indices[0] = 0; indices[1] = 3; indices[2] = 5; 49 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; 50 ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr); 51 row = 3; 52 indices[0] = 0; indices[1] = 1; indices[2] = 4; 53 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; 54 ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr); 55 row = 4; 56 indices[0] = 0; indices[1] = 1; indices[2] = 2; 57 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; 58 ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr); 59 break; 60 case 3: 61 row = 0; 62 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; 63 ierr = MatSetValues(A,1,&row,3,indices+1,values,INSERT_VALUES);CHKERRQ(ierr); 64 row = 1; 65 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; 66 ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr); 67 row = 2; 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 = 3; 71 values[0] = -2.0; values[1] = -2.0; values[2] = -2.0; values[3] = -1.0; 72 ierr = MatSetValues(A,1,&row,4,indices,values,INSERT_VALUES);CHKERRQ(ierr); 73 row = 4; 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 break; 77 78 default: 79 row = 0; 80 values[0] = -1.0; values[1] = 0.0; values[2] = 1.0; values[3] = 3.0; values[4] = 4.0; values[5] = -5.0; 81 ierr = MatSetValues(A,1,&row,N,indices,values,INSERT_VALUES);CHKERRQ(ierr); 82 row = 1; 83 ierr = MatSetValues(A,1,&row,3,indices,values,INSERT_VALUES);CHKERRQ(ierr); 84 row = 3; 85 ierr = MatSetValues(A,1,&row,1,indices+4,values+4,INSERT_VALUES);CHKERRQ(ierr); 86 row = 4; 87 ierr = MatSetValues(A,1,&row,2,indices+4,values+4,INSERT_VALUES);CHKERRQ(ierr); 88 } 89 } 90 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 91 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 92 ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 93 94 ierr = MatGetLocalSize(A, &m,&n);CHKERRQ(ierr); 95 ierr = VecCreate(PETSC_COMM_WORLD,&min);CHKERRQ(ierr); 96 ierr = VecSetSizes(min,m,PETSC_DECIDE);CHKERRQ(ierr); 97 ierr = VecSetFromOptions(min);CHKERRQ(ierr); 98 ierr = VecDuplicate(min,&max);CHKERRQ(ierr); 99 ierr = VecDuplicate(min,&maxabs);CHKERRQ(ierr); 100 ierr = VecDuplicate(min,&e);CHKERRQ(ierr); 101 102 /* MatGetRowMax() */ 103 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMax\n");CHKERRQ(ierr); 104 ierr = MatGetRowMax(A,max,NULL);CHKERRQ(ierr); 105 ierr = MatGetRowMax(A,max,imax);CHKERRQ(ierr); 106 ierr = VecView(max,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 107 ierr = VecGetLocalSize(max,&n);CHKERRQ(ierr); 108 ierr = PetscIntView(n,imax,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 109 110 /* MatGetRowMin() */ 111 ierr = MatScale(A,-1.0);CHKERRQ(ierr); 112 ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 113 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMin\n");CHKERRQ(ierr); 114 ierr = MatGetRowMin(A,min,NULL);CHKERRQ(ierr); 115 ierr = MatGetRowMin(A,min,imin);CHKERRQ(ierr); 116 117 ierr = VecWAXPY(e,1.0,max,min);CHKERRQ(ierr); /* e = max + min */ 118 ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr); 119 if (enorm > PETSC_MACHINE_EPSILON) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"max+min > PETSC_MACHINE_EPSILON "); 120 for (j = 0; j < n; j++) { 121 if (imin[j] != imax[j]) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"imin[%D] %D != imax %D",j,imin[j],imax[j]); 122 } 123 124 /* MatGetRowMaxAbs() */ 125 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMaxAbs\n");CHKERRQ(ierr); 126 ierr = MatGetRowMaxAbs(A,maxabs,NULL);CHKERRQ(ierr); 127 ierr = MatGetRowMaxAbs(A,maxabs,imaxabs);CHKERRQ(ierr); 128 ierr = VecView(maxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 129 ierr = PetscIntView(n,imaxabs,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 130 131 /* MatGetRowMinAbs() */ 132 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMinAbs\n");CHKERRQ(ierr); 133 ierr = MatGetRowMinAbs(A,min,NULL);CHKERRQ(ierr); 134 ierr = MatGetRowMinAbs(A,min,imin);CHKERRQ(ierr); 135 ierr = VecView(min,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 136 ierr = PetscIntView(n,imin,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 137 138 if (size == 1) { 139 /* Test MatGetRowMax, MatGetRowMin and MatGetRowMaxAbs for SeqDense and MPIBAIJ matrix */ 140 Mat Adense; 141 Vec max_d,maxabs_d; 142 ierr = VecDuplicate(min,&max_d);CHKERRQ(ierr); 143 ierr = VecDuplicate(min,&maxabs_d);CHKERRQ(ierr); 144 145 ierr = MatScale(A,-1.0);CHKERRQ(ierr); 146 ierr = MatConvert(A,MATDENSE,MAT_INITIAL_MATRIX,&Adense);CHKERRQ(ierr); 147 ierr = PetscPrintf(PETSC_COMM_WORLD,"MatGetRowMax for seqdense matrix\n");CHKERRQ(ierr); 148 ierr = MatGetRowMax(Adense,max_d,imax);CHKERRQ(ierr); 149 150 ierr = VecWAXPY(e,-1.0,max,max_d);CHKERRQ(ierr); /* e = -max + max_d */ 151 ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr); 152 if (enorm > PETSC_MACHINE_EPSILON) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"norm(-max + max_d) %g > PETSC_MACHINE_EPSILON",enorm); 153 154 ierr = MatScale(Adense,-1.0);CHKERRQ(ierr); 155 ierr = PetscPrintf(PETSC_COMM_WORLD,"MatGetRowMin for seqdense matrix\n");CHKERRQ(ierr); 156 ierr = MatGetRowMin(Adense,min,imin);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 SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"imin[%D] %D != imax %D for seqdense matrix",j,imin[j],imax[j]); 164 } 165 } 166 167 ierr = PetscPrintf(PETSC_COMM_WORLD,"MatGetRowMaxAbs for seqdense matrix\n");CHKERRQ(ierr); 168 ierr = MatGetRowMaxAbs(Adense,maxabs_d,imaxabs);CHKERRQ(ierr); 169 ierr = VecWAXPY(e,-1.0,maxabs,maxabs_d);CHKERRQ(ierr); /* e = -maxabs + maxabs_d */ 170 ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr); 171 if (enorm > PETSC_MACHINE_EPSILON) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"norm(-maxabs + maxabs_d) %g > PETSC_MACHINE_EPSILON",enorm); 172 173 ierr = MatDestroy(&Adense);CHKERRQ(ierr); 174 ierr = VecDestroy(&max_d);CHKERRQ(ierr); 175 ierr = VecDestroy(&maxabs_d);CHKERRQ(ierr); 176 } else { 177 /* BAIJ */ 178 Mat B; 179 Vec maxabsB; 180 PetscInt imaxabsB[M]; 181 ierr = MatConvert(A,MATMPIBAIJ,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr); 182 ierr = VecDuplicate(min,&maxabsB);CHKERRQ(ierr); 183 ierr = MatGetRowMaxAbs(A,maxabsB,NULL);CHKERRQ(ierr); 184 ierr = MatGetRowMaxAbs(A,maxabsB,imaxabsB);CHKERRQ(ierr); 185 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n MatGetRowMaxAbs for MPIBAIJ matrix\n");CHKERRQ(ierr); 186 ierr = VecWAXPY(e,-1.0,maxabs,maxabsB);CHKERRQ(ierr); /* e = -maxabs + maxabsB */ 187 ierr = VecNorm(e,NORM_INFINITY,&enorm);CHKERRQ(ierr); 188 if (enorm > PETSC_MACHINE_EPSILON) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"norm(-maxabs + maxabs_d) %g > PETSC_MACHINE_EPSILON",enorm); 189 190 for (j = 0; j < n; j++) { 191 if (imaxabs[j] != imaxabsB[j]) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"imaxabs[%D] %D != imaxabsB %D",j,imin[j],imax[j]); 192 } 193 ierr = MatDestroy(&B);CHKERRQ(ierr); 194 ierr = VecDestroy(&maxabsB);CHKERRQ(ierr); 195 } 196 197 ierr = VecDestroy(&min);CHKERRQ(ierr); 198 ierr = VecDestroy(&max);CHKERRQ(ierr); 199 ierr = VecDestroy(&maxabs);CHKERRQ(ierr); 200 ierr = VecDestroy(&e);CHKERRQ(ierr); 201 ierr = MatDestroy(&A);CHKERRQ(ierr); 202 ierr = PetscFinalize(); 203 return ierr; 204 } 205 206 /*TEST 207 208 test: 209 output_file: output/ex114.out 210 211 test: 212 suffix: 2 213 args: -testcase 1 214 output_file: output/ex114.out 215 216 test: 217 suffix: 3 218 args: -testcase 2 219 output_file: output/ex114_3.out 220 221 test: 222 suffix: 4 223 args: -testcase 3 224 output_file: output/ex114_4.out 225 226 test: 227 suffix: 5 228 nsize: 3 229 args: -testcase 0 230 output_file: output/ex114_5.out 231 232 test: 233 suffix: 6 234 nsize: 3 235 args: -testcase 1 236 output_file: output/ex114_6.out 237 238 test: 239 suffix: 7 240 nsize: 3 241 args: -testcase 2 242 output_file: output/ex114_7.out 243 244 test: 245 suffix: 8 246 nsize: 3 247 args: -testcase 3 248 output_file: output/ex114_8.out 249 250 TEST*/ 251