1 2 static char help[] = "Test Hypre matrix APIs\n"; 3 4 #include <petscmathypre.h> 5 6 int main(int argc,char **args) 7 { 8 Mat A, B, C; 9 PetscReal err; 10 PetscInt i,j,M = 20; 11 PetscMPIInt NP; 12 MPI_Comm comm; 13 PetscInt *rows; 14 15 CHKERRQ(PetscInitialize(&argc,&args,(char*)0,help)); 16 comm = PETSC_COMM_WORLD; 17 CHKERRMPI(MPI_Comm_size(comm,&NP)); 18 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 19 PetscCheckFalse(M < 6,PETSC_COMM_WORLD,PETSC_ERR_SUP,"Matrix has to have more than 6 columns"); 20 /* Hypre matrix */ 21 CHKERRQ(MatCreate(comm,&B)); 22 CHKERRQ(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,M,M)); 23 CHKERRQ(MatSetType(B,MATHYPRE)); 24 CHKERRQ(MatHYPRESetPreallocation(B,9,NULL,9,NULL)); 25 26 /* PETSc AIJ matrix */ 27 CHKERRQ(MatCreate(comm,&A)); 28 CHKERRQ(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,M)); 29 CHKERRQ(MatSetType(A,MATAIJ)); 30 CHKERRQ(MatSeqAIJSetPreallocation(A,9,NULL)); 31 CHKERRQ(MatMPIAIJSetPreallocation(A,9,NULL,9,NULL)); 32 33 /*Set Values */ 34 for (i=0; i<M; i++) { 35 PetscInt cols[] = {0,1,2,3,4,5}; 36 PetscScalar vals[6] = {0}; 37 PetscScalar value[] = {100}; 38 for (j=0; j<6; j++) 39 vals[j] = ((PetscReal)j)/NP; 40 41 CHKERRQ(MatSetValues(B,1,&i,6,cols,vals,ADD_VALUES)); 42 CHKERRQ(MatSetValues(B,1,&i,1,&i,value,ADD_VALUES)); 43 CHKERRQ(MatSetValues(A,1,&i,6,cols,vals,ADD_VALUES)); 44 CHKERRQ(MatSetValues(A,1,&i,1,&i,value,ADD_VALUES)); 45 } 46 47 /* MAT_FLUSH_ASSEMBLY currently not supported */ 48 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 49 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 50 CHKERRQ(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); 51 CHKERRQ(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY)); 52 53 /* Compare A and B */ 54 CHKERRQ(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 55 CHKERRQ(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 56 CHKERRQ(MatNorm(C,NORM_INFINITY,&err)); 57 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)B),PETSC_ERR_PLIB,"Error MatSetValues %g",err); 58 CHKERRQ(MatDestroy(&C)); 59 60 /* MatZeroRows */ 61 CHKERRQ(PetscMalloc1(M, &rows)); 62 for (i=0; i<M; i++) rows[i] = i; 63 CHKERRQ(MatZeroRows(B, M, rows, 10.0, NULL, NULL)); 64 CHKERRQ(MatSetOption(A, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE)); 65 CHKERRQ(MatZeroRows(A, M, rows, 10.0,NULL, NULL)); 66 CHKERRQ(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 67 CHKERRQ(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 68 CHKERRQ(MatNorm(C,NORM_INFINITY,&err)); 69 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)B),PETSC_ERR_PLIB,"Error MatZeroRows %g",err); 70 CHKERRQ(MatDestroy(&C)); 71 CHKERRQ(PetscFree(rows)); 72 73 /* Test MatZeroEntries */ 74 CHKERRQ(MatZeroEntries(B)); 75 CHKERRQ(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 76 CHKERRQ(MatNorm(C,NORM_INFINITY,&err)); 77 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatZeroEntries %g",err); 78 CHKERRQ(MatDestroy(&C)); 79 80 /* Insert Values */ 81 for (i=0; i<M; i++) { 82 PetscInt cols[] = {0,1,2,3,4,5}; 83 PetscScalar vals[6] = {0}; 84 PetscScalar value[] = {100}; 85 86 for (j=0; j<6; j++) 87 vals[j] = ((PetscReal)j)/NP; 88 89 CHKERRQ(MatSetValues(B,1,&i,6,cols,vals,INSERT_VALUES)); 90 CHKERRQ(MatSetValues(B,1,&i,1,&i,value,INSERT_VALUES)); 91 CHKERRQ(MatSetValues(A,1,&i,6,cols,vals,INSERT_VALUES)); 92 CHKERRQ(MatSetValues(A,1,&i,1,&i,value,INSERT_VALUES)); 93 } 94 95 /* MAT_FLUSH_ASSEMBLY currently not supported */ 96 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 97 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 98 CHKERRQ(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); 99 CHKERRQ(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY)); 100 101 /* Rows are not sorted with HYPRE so we need an intermediate sort 102 They use a temporary buffer, so we can sort inplace the const memory */ 103 { 104 const PetscInt *idxA,*idxB; 105 const PetscScalar *vA, *vB; 106 PetscInt rstart, rend, nzA, nzB; 107 PetscInt cols[] = {0,1,2,3,4,-5}; 108 PetscInt *rows; 109 PetscScalar *valuesA, *valuesB; 110 PetscBool flg; 111 112 CHKERRQ(MatGetOwnershipRange(A,&rstart,&rend)); 113 for (i=rstart; i<rend; i++) { 114 CHKERRQ(MatGetRow(A,i,&nzA,&idxA,&vA)); 115 CHKERRQ(MatGetRow(B,i,&nzB,&idxB,&vB)); 116 PetscCheckFalse(nzA!=nzB,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT, nzA-nzB); 117 CHKERRQ(PetscSortIntWithScalarArray(nzB,(PetscInt*)idxB,(PetscScalar*)vB)); 118 CHKERRQ(PetscArraycmp(idxA,idxB,nzA,&flg)); 119 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT " (indices)",i); 120 CHKERRQ(PetscArraycmp(vA,vB,nzA,&flg)); 121 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT " (values)",i); 122 CHKERRQ(MatRestoreRow(A,i,&nzA,&idxA,&vA)); 123 CHKERRQ(MatRestoreRow(B,i,&nzB,&idxB,&vB)); 124 } 125 126 CHKERRQ(MatGetOwnershipRange(A,&rstart,&rend)); 127 CHKERRQ(PetscCalloc3((rend-rstart)*6,&valuesA,(rend-rstart)*6,&valuesB,rend-rstart,&rows)); 128 for (i=rstart; i<rend; i++) rows[i-rstart] =i; 129 130 CHKERRQ(MatGetValues(A,rend-rstart,rows,6,cols,valuesA)); 131 CHKERRQ(MatGetValues(B,rend-rstart,rows,6,cols,valuesB)); 132 133 for (i=0; i<(rend-rstart); i++) { 134 CHKERRQ(PetscArraycmp(valuesA + 6*i,valuesB + 6*i,6,&flg)); 135 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetValues %" PetscInt_FMT,i + rstart); 136 } 137 CHKERRQ(PetscFree3(valuesA,valuesB,rows)); 138 } 139 140 /* Compare A and B */ 141 CHKERRQ(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 142 CHKERRQ(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 143 CHKERRQ(MatNorm(C,NORM_INFINITY,&err)); 144 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)B),PETSC_ERR_PLIB,"Error MatSetValues with INSERT_VALUES %g",err); 145 146 CHKERRQ(MatDestroy(&A)); 147 CHKERRQ(MatDestroy(&B)); 148 CHKERRQ(MatDestroy(&C)); 149 150 CHKERRQ(PetscFinalize()); 151 return 0; 152 } 153 154 /*TEST 155 156 build: 157 requires: hypre 158 159 test: 160 suffix: 1 161 requires: !defined(PETSC_HAVE_HYPRE_DEVICE) 162 163 test: 164 suffix: 2 165 requires: !defined(PETSC_HAVE_HYPRE_DEVICE) 166 output_file: output/ex225_1.out 167 nsize: 2 168 169 TEST*/ 170