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 PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 16 comm = PETSC_COMM_WORLD; 17 PetscCallMPI(MPI_Comm_size(comm,&NP)); 18 PetscCall(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 PetscCall(MatCreate(comm,&B)); 22 PetscCall(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,M,M)); 23 PetscCall(MatSetType(B,MATHYPRE)); 24 PetscCall(MatHYPRESetPreallocation(B,9,NULL,9,NULL)); 25 26 /* PETSc AIJ matrix */ 27 PetscCall(MatCreate(comm,&A)); 28 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,M)); 29 PetscCall(MatSetType(A,MATAIJ)); 30 PetscCall(MatSeqAIJSetPreallocation(A,9,NULL)); 31 PetscCall(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 PetscCall(MatSetValues(B,1,&i,6,cols,vals,ADD_VALUES)); 42 PetscCall(MatSetValues(B,1,&i,1,&i,value,ADD_VALUES)); 43 PetscCall(MatSetValues(A,1,&i,6,cols,vals,ADD_VALUES)); 44 PetscCall(MatSetValues(A,1,&i,1,&i,value,ADD_VALUES)); 45 } 46 47 /* MAT_FLUSH_ASSEMBLY currently not supported */ 48 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 49 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 50 PetscCall(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); 51 PetscCall(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY)); 52 53 /* Compare A and B */ 54 PetscCall(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 55 PetscCall(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 56 PetscCall(MatNorm(C,NORM_INFINITY,&err)); 57 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)B),PETSC_ERR_PLIB,"Error MatSetValues %g",err); 58 PetscCall(MatDestroy(&C)); 59 60 /* MatZeroRows */ 61 PetscCall(PetscMalloc1(M, &rows)); 62 for (i=0; i<M; i++) rows[i] = i; 63 PetscCall(MatZeroRows(B, M, rows, 10.0, NULL, NULL)); 64 PetscCall(MatSetOption(A, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE)); 65 PetscCall(MatZeroRows(A, M, rows, 10.0,NULL, NULL)); 66 PetscCall(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 67 PetscCall(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 68 PetscCall(MatNorm(C,NORM_INFINITY,&err)); 69 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)B),PETSC_ERR_PLIB,"Error MatZeroRows %g",err); 70 PetscCall(MatDestroy(&C)); 71 PetscCall(PetscFree(rows)); 72 73 /* Test MatZeroEntries */ 74 PetscCall(MatZeroEntries(B)); 75 PetscCall(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 76 PetscCall(MatNorm(C,NORM_INFINITY,&err)); 77 PetscCheckFalse(err > PETSC_SMALL,PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatZeroEntries %g",err); 78 PetscCall(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 PetscCall(MatSetValues(B,1,&i,6,cols,vals,INSERT_VALUES)); 90 PetscCall(MatSetValues(B,1,&i,1,&i,value,INSERT_VALUES)); 91 PetscCall(MatSetValues(A,1,&i,6,cols,vals,INSERT_VALUES)); 92 PetscCall(MatSetValues(A,1,&i,1,&i,value,INSERT_VALUES)); 93 } 94 95 /* MAT_FLUSH_ASSEMBLY currently not supported */ 96 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 97 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 98 PetscCall(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); 99 PetscCall(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 PetscCall(MatGetOwnershipRange(A,&rstart,&rend)); 113 for (i=rstart; i<rend; i++) { 114 PetscCall(MatGetRow(A,i,&nzA,&idxA,&vA)); 115 PetscCall(MatGetRow(B,i,&nzB,&idxB,&vB)); 116 PetscCheckFalse(nzA!=nzB,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT, nzA-nzB); 117 PetscCall(PetscSortIntWithScalarArray(nzB,(PetscInt*)idxB,(PetscScalar*)vB)); 118 PetscCall(PetscArraycmp(idxA,idxB,nzA,&flg)); 119 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT " (indices)",i); 120 PetscCall(PetscArraycmp(vA,vB,nzA,&flg)); 121 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error MatGetRow %" PetscInt_FMT " (values)",i); 122 PetscCall(MatRestoreRow(A,i,&nzA,&idxA,&vA)); 123 PetscCall(MatRestoreRow(B,i,&nzB,&idxB,&vB)); 124 } 125 126 PetscCall(MatGetOwnershipRange(A,&rstart,&rend)); 127 PetscCall(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 PetscCall(MatGetValues(A,rend-rstart,rows,6,cols,valuesA)); 131 PetscCall(MatGetValues(B,rend-rstart,rows,6,cols,valuesB)); 132 133 for (i=0; i<(rend-rstart); i++) { 134 PetscCall(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 PetscCall(PetscFree3(valuesA,valuesB,rows)); 138 } 139 140 /* Compare A and B */ 141 PetscCall(MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C)); 142 PetscCall(MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN)); 143 PetscCall(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 PetscCall(MatDestroy(&A)); 147 PetscCall(MatDestroy(&B)); 148 PetscCall(MatDestroy(&C)); 149 150 PetscCall(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