Actual source code: seqhashmatsetvalues.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 and bs > 1
  6:             TYPE_SBAIJ for SBAIJ

  8: */
  9: static PetscErrorCode PetscConcat(MatSetValues_Seq_Hash, TYPE_BS)(Mat A, PetscInt m, const PetscInt *rows, PetscInt n, const PetscInt *cols, const PetscScalar *values, InsertMode addv)
 10: {
 11:   PetscConcat(Mat_Seq, TYPE) *a = (PetscConcat(Mat_Seq, TYPE) *)A->data;
 12: #if defined(TYPE_BS_ON)
 13:   PetscInt bs;
 14: #endif

 16:   PetscFunctionBegin;
 17: #if defined(TYPE_BS_ON)
 18:   PetscCall(MatGetBlockSize(A, &bs));
 19: #endif
 20:   for (PetscInt r = 0; r < m; ++r) {
 21:     PetscHashIJKey key;
 22:     PetscBool      missing;
 23:     PetscScalar    value;
 24: #if defined(TYPE_BS_ON)
 25:     PetscHashIJKey bkey;
 26: #endif

 28:     key.i = rows[r];
 29: #if defined(TYPE_BS_ON)
 30:     bkey.i = key.i / bs;
 31: #endif
 32:     if (key.i < 0) continue;
 33:     for (PetscInt c = 0; c < n; ++c) {
 34:       key.j = cols[c];
 35: #if defined(TYPE_BS_ON)
 36:       bkey.j = key.j / bs;
 37:   #if defined(TYPE_SBAIJ)
 38:       if (bkey.j < bkey.i) continue;
 39:   #else
 40:       if (key.j < 0) continue;
 41:   #endif
 42: #else
 43:   #if defined(TYPE_SBAIJ)
 44:       if (key.j < key.i) continue;
 45:   #else
 46:       if (key.j < 0) continue;
 47:   #endif
 48: #endif
 49:       value = values ? ((a->roworiented) ? values[r * n + c] : values[r + m * c]) : 0;
 50:       switch (addv) {
 51:       case INSERT_VALUES:
 52:         PetscCall(PetscHMapIJVQuerySet(a->ht, key, value, &missing));
 53:         break;
 54:       case ADD_VALUES:
 55:         PetscCall(PetscHMapIJVQueryAdd(a->ht, key, value, &missing));
 56:         break;
 57:       default:
 58:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "InsertMode not supported");
 59:       }
 60:       if (missing) ++a->dnz[key.i];
 61: #if defined(TYPE_BS_ON)
 62:       PetscCall(PetscHSetIJQueryAdd(a->bht, bkey, &missing));
 63:       if (missing) ++a->bdnz[bkey.i];
 64: #endif
 65:     }
 66:   }
 67:   PetscFunctionReturn(PETSC_SUCCESS);
 68: }