MatCopyHashToXAIJ_Seq_Hash(Mat A,Mat B)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 */
MatAssemblyEnd_Seq_Hash(Mat A,MatAssemblyType type)70 static PetscErrorCode MatAssemblyEnd_Seq_Hash(Mat A, MatAssemblyType type)
71 {
72 PetscFunctionBegin;
73 PetscCall(MatCopyHashToXAIJ(A, A));
74 PetscFunctionReturn(PETSC_SUCCESS);
75 }
76
MatDestroy_Seq_Hash(Mat A)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
MatZeroEntries_Seq_Hash(Mat A)98 static PetscErrorCode MatZeroEntries_Seq_Hash(Mat A)
99 {
100 PetscFunctionBegin;
101 PetscFunctionReturn(PETSC_SUCCESS);
102 }
103
MatSetRandom_Seq_Hash(Mat A,PetscRandom r)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
MatSetUp_Seq_Hash(Mat A)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