Actual source code: seqhashmat.h

  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

 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 */

 20:     A->ops[0]      = a->cops;
 21:     A->hash_active = PETSC_FALSE;
 22:   }

 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, B->structure_only ? 0 : n, &values));
 36:   rowstarts[0] = 0;
 37:   for (PetscInt i = 0; i < m; i++) rowstarts[i + 1] = rowstarts[i] + a->dnz[i];

 39:   PetscHashIterBegin(a->ht, hi);
 40:   while (!PetscHashIterAtEnd(a->ht, hi)) {
 41:     PetscHashIterGetKey(a->ht, hi, key);
 42:     cols[rowstarts[key.i]] = key.j;
 43:     if (!B->structure_only) {
 44:       PetscHashIterGetVal(a->ht, hi, value);
 45:       values[rowstarts[key.i]] = value;
 46:     }
 47:     rowstarts[key.i]++;
 48:     PetscHashIterNext(a->ht, hi);
 49:   }
 50:   if (A == B) PetscCall(PetscHMapIJVDestroy(&a->ht));

 52:   for (PetscInt i = 0, start = 0; i < m; i++) {
 53:     PetscCall(MatSetValues(B, 1, &i, a->dnz[i], PetscSafePointerPlusOffset(cols, start), PetscSafePointerPlusOffset(values, start), B->insertmode));
 54:     start += a->dnz[i];
 55:   }
 56:   PetscCall(PetscFree3(cols, rowstarts, values));
 57:   if (A == B) PetscCall(PetscFree(a->dnz));
 58: #if defined(TYPE_BS_ON)
 59:   if (bs > 1 && A == B) PetscCall(PetscFree(a->bdnz));
 60: #endif
 61:   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
 62:   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
 63:   PetscFunctionReturn(PETSC_SUCCESS);
 64: }

 66: /*
 67:    used by SEQAIJ, BAIJ and SBAIJ to reduce code duplication

 69:      define TYPE to AIJ BAIJ or SBAIJ
 70:             TYPE_BS_ON for BAIJ and SBAIJ

 72: */
 73: static PetscErrorCode MatAssemblyEnd_Seq_Hash(Mat A, MatAssemblyType type)
 74: {
 75:   PetscFunctionBegin;
 76:   PetscCall(MatCopyHashToXAIJ(A, A));
 77:   PetscFunctionReturn(PETSC_SUCCESS);
 78: }

 80: static PetscErrorCode MatDestroy_Seq_Hash(Mat A)
 81: {
 82:   PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data;
 83: #if defined(TYPE_BS_ON)
 84:   PetscInt bs;
 85: #endif

 87:   PetscFunctionBegin;
 88:   PetscCall(PetscHMapIJVDestroy(&a->ht));
 89:   PetscCall(PetscFree(a->dnz));
 90: #if defined(TYPE_BS_ON)
 91:   PetscCall(MatGetBlockSize(A, &bs));
 92:   if (bs > 1) {
 93:     PetscCall(PetscFree(a->bdnz));
 94:     PetscCall(PetscHSetIJDestroy(&a->bht));
 95:   }
 96: #endif
 97:   PetscCall((*a->cops.destroy)(A));
 98:   PetscFunctionReturn(PETSC_SUCCESS);
 99: }

101: static PetscErrorCode MatZeroEntries_Seq_Hash(Mat A)
102: {
103:   PetscFunctionBegin;
104:   PetscFunctionReturn(PETSC_SUCCESS);
105: }

107: static PetscErrorCode MatSetRandom_Seq_Hash(Mat A, PetscRandom r)
108: {
109:   PetscFunctionBegin;
110:   SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Must set preallocation first");
111:   PetscFunctionReturn(PETSC_SUCCESS);
112: }

114: static PetscErrorCode MatSetUp_Seq_Hash(Mat A)
115: {
116:   PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data;
117:   PetscInt m;
118: #if defined(TYPE_BS_ON)
119:   PetscInt bs;
120: #endif

122:   PetscFunctionBegin;
123:   PetscCall(PetscInfo(A, "Using hash-based MatSetValues() for MATSEQ" PetscStringize(TYPE) " because no preallocation provided\n"));
124:   PetscCall(PetscLayoutSetUp(A->rmap));
125:   PetscCall(PetscLayoutSetUp(A->cmap));
126:   if (A->rmap->bs < 1) A->rmap->bs = 1;
127:   if (A->cmap->bs < 1) A->cmap->bs = 1;

129:   PetscCall(MatGetLocalSize(A, &m, NULL));
130:   PetscCall(PetscHMapIJVCreate(&a->ht));
131:   PetscCall(PetscCalloc1(m, &a->dnz));
132: #if defined(TYPE_BS_ON)
133:   PetscCall(MatGetBlockSize(A, &bs));
134:   if (bs > 1) {
135:     PetscCall(PetscCalloc1(m / bs, &a->bdnz));
136:     PetscCall(PetscHSetIJCreate(&a->bht));
137:   }
138: #endif

140:   /* keep a record of the operations so they can be reset when the hash handling is complete */
141:   a->cops                = A->ops[0];
142:   A->ops->assemblybegin  = NULL;
143:   A->ops->assemblyend    = MatAssemblyEnd_Seq_Hash;
144:   A->ops->destroy        = MatDestroy_Seq_Hash;
145:   A->ops->zeroentries    = MatZeroEntries_Seq_Hash;
146:   A->ops->setrandom      = MatSetRandom_Seq_Hash;
147:   A->ops->copyhashtoxaij = MatCopyHashToXAIJ_Seq_Hash;
148: #if defined(TYPE_BS_ON)
149:   if (bs > 1) A->ops->setvalues = MatSetValues_Seq_Hash_BS;
150:   else
151: #endif
152:     A->ops->setvalues = MatSetValues_Seq_Hash;
153:   A->ops->setvaluesblocked = NULL;

155:   A->preallocated = PETSC_TRUE;
156:   A->hash_active  = PETSC_TRUE;
157:   PetscFunctionReturn(PETSC_SUCCESS);
158: }