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: PetscFunctionBegin;
19: *nn = n;
20: if (!ia) PetscFunctionReturn(PETSC_SUCCESS);
22: PetscCall(PetscCalloc1(n + 1, &collengths));
23: PetscCall(PetscMalloc1(n + 1, &cia));
24: PetscCall(PetscMalloc1(a->nz + 1, &cja));
25: PetscCall(PetscMalloc1(a->nz + 1, &cspidx));
27: totalslices = A->rmap->n / 8 + ((A->rmap->n & 0x07) ? 1 : 0); /* floor(n/8) */
28: for (i = 0; i < totalslices; i++) { /* loop over slices */
29: for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
30: isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
31: if (isnonzero) collengths[a->colidx[j]]++;
32: }
33: }
35: cia[0] = oshift;
36: for (i = 0; i < n; i++) cia[i + 1] = cia[i] + collengths[i];
37: PetscCall(PetscArrayzero(collengths, n));
39: for (i = 0; i < totalslices; i++) { /* loop over slices */
40: for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
41: isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
42: if (isnonzero) {
43: col = a->colidx[j];
44: cspidx[cia[col] + collengths[col] - oshift] = j; /* index of a->colidx */
45: cja[cia[col] + collengths[col] - oshift] = 8 * i + row + oshift; /* row index */
46: collengths[col]++;
47: }
48: }
49: }
51: PetscCall(PetscFree(collengths));
52: *ia = cia;
53: *ja = cja;
54: *spidx = cspidx;
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }
58: PetscErrorCode MatRestoreColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *n, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done)
59: {
60: PetscFunctionBegin;
62: if (!ia) PetscFunctionReturn(PETSC_SUCCESS);
63: PetscCall(PetscFree(*ia));
64: PetscCall(PetscFree(*ja));
65: PetscCall(PetscFree(*spidx));
66: PetscFunctionReturn(PETSC_SUCCESS);
67: }