Actual source code: fdsell.c

  1: #include <../src/mat/impls/sell/seq/sell.h>
  2: #include <../src/mat/impls/aij/seq/aij.h>
  3: #include <petsc/private/isimpl.h>

  5: /*
  6:  MatGetColumnIJ_SeqSELL_Color() and MatRestoreColumnIJ_SeqSELL_Color() are customized from
  7:  MatGetColumnIJ_SeqSELL() and MatRestoreColumnIJ_SeqSELL() by adding an output
  8:  spidx[], index of a->a, to be used in MatTransposeColoringCreate_SeqSELL() and MatFDColoringCreate_SeqSELL()
  9: */
 10: PetscErrorCode MatGetColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *nn, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done)
 11: {
 12:   Mat_SeqSELL *a = (Mat_SeqSELL *)A->data;
 13:   PetscInt     i, j, *collengths, *cia, *cja, n = A->cmap->n, totalslices;
 14:   PetscInt     row, col;
 15:   PetscInt    *cspidx;
 16:   PetscBool    isnonzero;

 18:   *nn = n;
 19:   if (!ia) return 0;

 21:   PetscCalloc1(n + 1, &collengths);
 22:   PetscMalloc1(n + 1, &cia);
 23:   PetscMalloc1(a->nz + 1, &cja);
 24:   PetscMalloc1(a->nz + 1, &cspidx);

 26:   totalslices = A->rmap->n / 8 + ((A->rmap->n & 0x07) ? 1 : 0); /* floor(n/8) */
 27:   for (i = 0; i < totalslices; i++) {                           /* loop over slices */
 28:     for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
 29:       isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
 30:       if (isnonzero) collengths[a->colidx[j]]++;
 31:     }
 32:   }

 34:   cia[0] = oshift;
 35:   for (i = 0; i < n; i++) cia[i + 1] = cia[i] + collengths[i];
 36:   PetscArrayzero(collengths, n);

 38:   for (i = 0; i < totalslices; i++) { /* loop over slices */
 39:     for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
 40:       isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
 41:       if (isnonzero) {
 42:         col                                         = a->colidx[j];
 43:         cspidx[cia[col] + collengths[col] - oshift] = j;                    /* index of a->colidx */
 44:         cja[cia[col] + collengths[col] - oshift]    = 8 * i + row + oshift; /* row index */
 45:         collengths[col]++;
 46:       }
 47:     }
 48:   }

 50:   PetscFree(collengths);
 51:   *ia    = cia;
 52:   *ja    = cja;
 53:   *spidx = cspidx;
 54:   return 0;
 55: }

 57: PetscErrorCode MatRestoreColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *n, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done)
 58: {

 60:   if (!ia) return 0;
 61:   PetscFree(*ia);
 62:   PetscFree(*ja);
 63:   PetscFree(*spidx);
 64:   return 0;
 65: }