1 /* 2 used by SEQAIJ, BAIJ and SBAIJ to reduce code duplication 3 4 define TYPE to AIJ BAIJ or SBAIJ 5 TYPE_BS_ON for BAIJ and SBAIJ 6 7 */ 8 static PetscErrorCode MatAssemblyEnd_Seq_Hash(Mat A, MatAssemblyType type) 9 { 10 PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data; 11 PetscHashIter hi; 12 PetscHashIJKey key; 13 PetscScalar value, *values; 14 PetscInt m, n, *cols, *rowstarts; 15 #if defined(TYPE_BS_ON) 16 PetscInt bs; 17 #endif 18 19 PetscFunctionBegin; 20 #if defined(TYPE_BS_ON) 21 PetscCall(MatGetBlockSize(A, &bs)); 22 if (bs > 1) PetscCall(PetscHSetIJDestroy(&a->bht)); 23 #endif 24 A->preallocated = PETSC_FALSE; /* this was set to true for the MatSetValues_Hash() to work */ 25 26 A->ops[0] = a->cops; 27 A->hash_active = PETSC_FALSE; 28 29 /* move values from hash format to matrix type format */ 30 PetscCall(MatGetSize(A, &m, NULL)); 31 #if defined(TYPE_BS_ON) 32 if (bs > 1) PetscCall(PetscConcat(PetscConcat(MatSeq, TYPE), SetPreallocation)(A, bs, PETSC_DETERMINE, a->bdnz)); 33 else PetscCall(PetscConcat(PetscConcat(MatSeq, TYPE), SetPreallocation)(A, 1, PETSC_DETERMINE, a->dnz)); 34 #else 35 PetscCall(MatSeqAIJSetPreallocation(A, PETSC_DETERMINE, a->dnz)); 36 #endif 37 PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE)); 38 PetscCall(PetscHMapIJVGetSize(a->ht, &n)); 39 /* do not need PetscShmgetAllocateArray() since arrays are temporary */ 40 PetscCall(PetscMalloc3(n, &cols, m + 1, &rowstarts, n, &values)); 41 rowstarts[0] = 0; 42 for (PetscInt i = 0; i < m; i++) rowstarts[i + 1] = rowstarts[i] + a->dnz[i]; 43 44 PetscHashIterBegin(a->ht, hi); 45 while (!PetscHashIterAtEnd(a->ht, hi)) { 46 PetscHashIterGetKey(a->ht, hi, key); 47 PetscHashIterGetVal(a->ht, hi, value); 48 cols[rowstarts[key.i]] = key.j; 49 values[rowstarts[key.i]++] = value; 50 PetscHashIterNext(a->ht, hi); 51 } 52 PetscCall(PetscHMapIJVDestroy(&a->ht)); 53 54 for (PetscInt i = 0, start = 0; i < m; i++) { 55 PetscCall(MatSetValues(A, 1, &i, a->dnz[i], PetscSafePointerPlusOffset(cols, start), PetscSafePointerPlusOffset(values, start), A->insertmode)); 56 start += a->dnz[i]; 57 } 58 PetscCall(PetscFree3(cols, rowstarts, values)); 59 PetscCall(PetscFree(a->dnz)); 60 #if defined(TYPE_BS_ON) 61 if (bs > 1) PetscCall(PetscFree(a->bdnz)); 62 #endif 63 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 64 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 65 PetscFunctionReturn(PETSC_SUCCESS); 66 } 67 68 static PetscErrorCode MatDestroy_Seq_Hash(Mat A) 69 { 70 PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data; 71 #if defined(TYPE_BS_ON) 72 PetscInt bs; 73 #endif 74 75 PetscFunctionBegin; 76 PetscCall(PetscHMapIJVDestroy(&a->ht)); 77 PetscCall(PetscFree(a->dnz)); 78 #if defined(TYPE_BS_ON) 79 PetscCall(MatGetBlockSize(A, &bs)); 80 if (bs > 1) { 81 PetscCall(PetscFree(a->bdnz)); 82 PetscCall(PetscHSetIJDestroy(&a->bht)); 83 } 84 #endif 85 PetscCall((*a->cops.destroy)(A)); 86 PetscFunctionReturn(PETSC_SUCCESS); 87 } 88 89 static PetscErrorCode MatZeroEntries_Seq_Hash(Mat A) 90 { 91 PetscFunctionBegin; 92 PetscFunctionReturn(PETSC_SUCCESS); 93 } 94 95 static PetscErrorCode MatSetRandom_Seq_Hash(Mat A, PetscRandom r) 96 { 97 PetscFunctionBegin; 98 SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Must set preallocation first"); 99 PetscFunctionReturn(PETSC_SUCCESS); 100 } 101 102 static PetscErrorCode MatSetUp_Seq_Hash(Mat A) 103 { 104 PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data; 105 PetscInt m; 106 #if defined(TYPE_BS_ON) 107 PetscInt bs; 108 #endif 109 110 PetscFunctionBegin; 111 PetscCall(PetscInfo(A, "Using hash-based MatSetValues() for MATSEQAIJ because no preallocation provided\n")); 112 PetscCall(PetscLayoutSetUp(A->rmap)); 113 PetscCall(PetscLayoutSetUp(A->cmap)); 114 if (A->rmap->bs < 1) A->rmap->bs = 1; 115 if (A->cmap->bs < 1) A->cmap->bs = 1; 116 117 PetscCall(MatGetLocalSize(A, &m, NULL)); 118 PetscCall(PetscHMapIJVCreate(&a->ht)); 119 PetscCall(PetscCalloc1(m, &a->dnz)); 120 #if defined(TYPE_BS_ON) 121 PetscCall(MatGetBlockSize(A, &bs)); 122 if (bs > 1) { 123 PetscCall(PetscCalloc1(m / bs, &a->bdnz)); 124 PetscCall(PetscHSetIJCreate(&a->bht)); 125 } 126 #endif 127 128 /* keep a record of the operations so they can be reset when the hash handling is complete */ 129 a->cops = A->ops[0]; 130 A->ops->assemblybegin = NULL; 131 A->ops->assemblyend = MatAssemblyEnd_Seq_Hash; 132 A->ops->destroy = MatDestroy_Seq_Hash; 133 A->ops->zeroentries = MatZeroEntries_Seq_Hash; 134 A->ops->setrandom = MatSetRandom_Seq_Hash; 135 #if defined(TYPE_BS_ON) 136 if (bs > 1) A->ops->setvalues = MatSetValues_Seq_Hash_BS; 137 else 138 #endif 139 A->ops->setvalues = MatSetValues_Seq_Hash; 140 A->ops->setvaluesblocked = NULL; 141 142 A->preallocated = PETSC_TRUE; 143 A->hash_active = PETSC_TRUE; 144 PetscFunctionReturn(PETSC_SUCCESS); 145 } 146