Actual source code: seqhashmat.h

  1: /*
  2:    used by SEQAIJ, BAIJ and SBAIJ to reduce code duplication

  4:      define TYPE to AIJ BAIJ or SBAIJ
  5:             TYPE_BS_ON for BAIJ and SBAIJ

  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

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

 26:   PetscCall(PetscMemcpy(&A->ops, &a->cops, sizeof(*(A->ops))));
 27:   A->hash_active = PETSC_FALSE;

 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:   PetscCall(PetscMalloc3(n, &cols, m + 1, &rowstarts, n, &values));
 40:   rowstarts[0] = 0;
 41:   for (PetscInt i = 0; i < m; i++) rowstarts[i + 1] = rowstarts[i] + a->dnz[i];

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

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

 67: static PetscErrorCode MatDestroy_Seq_Hash(Mat A)
 68: {
 69:   PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data;
 70: #if defined(TYPE_BS_ON)
 71:   PetscInt bs;
 72: #endif

 74:   PetscFunctionBegin;
 75:   PetscCall(PetscHMapIJVDestroy(&a->ht));
 76:   PetscCall(PetscFree(a->dnz));
 77: #if defined(TYPE_BS_ON)
 78:   PetscCall(MatGetBlockSize(A, &bs));
 79:   if (bs > 1) {
 80:     PetscCall(PetscFree(a->bdnz));
 81:     PetscCall(PetscHSetIJDestroy(&a->bht));
 82:   }
 83: #endif
 84:   PetscCall((*a->cops.destroy)(A));
 85:   PetscFunctionReturn(PETSC_SUCCESS);
 86: }

 88: static PetscErrorCode MatZeroEntries_Seq_Hash(Mat A)
 89: {
 90:   PetscFunctionBegin;
 91:   PetscFunctionReturn(PETSC_SUCCESS);
 92: }

 94: static PetscErrorCode MatSetRandom_Seq_Hash(Mat A, PetscRandom r)
 95: {
 96:   PetscFunctionBegin;
 97:   SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Must set preallocation first");
 98:   PetscFunctionReturn(PETSC_SUCCESS);
 99: }

101: static PetscErrorCode MatSetUp_Seq_Hash(Mat A)
102: {
103:   PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data;
104:   PetscInt m;
105: #if defined(TYPE_BS_ON)
106:   PetscInt bs;
107: #endif

109:   PetscFunctionBegin;
110:   PetscCall(PetscInfo(A, "Using hash-based MatSetValues() for MATSEQAIJ because no preallocation provided\n"));
111:   PetscCall(PetscLayoutSetUp(A->rmap));
112:   PetscCall(PetscLayoutSetUp(A->cmap));
113:   if (A->rmap->bs < 1) A->rmap->bs = 1;
114:   if (A->cmap->bs < 1) A->cmap->bs = 1;

116:   PetscCall(MatGetLocalSize(A, &m, NULL));
117:   PetscCall(PetscHMapIJVCreate(&a->ht));
118:   PetscCall(PetscCalloc1(m, &a->dnz));
119: #if defined(TYPE_BS_ON)
120:   PetscCall(MatGetBlockSize(A, &bs));
121:   if (bs > 1) {
122:     PetscCall(PetscCalloc1(m / bs, &a->bdnz));
123:     PetscCall(PetscHSetIJCreate(&a->bht));
124:   }
125: #endif

127:   /* keep a record of the operations so they can be reset when the hash handling is complete */
128:   PetscCall(PetscMemcpy(&a->cops, &A->ops, sizeof(*(A->ops))));

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;

142:   A->preallocated = PETSC_TRUE;
143:   A->hash_active  = PETSC_TRUE;
144:   PetscFunctionReturn(PETSC_SUCCESS);
145: }