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