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