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:   A->ops[0]      = a->cops;
     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:   /* 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];

     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));

     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: }

     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

     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: }

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

     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: }

    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

    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;

    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

    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;

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