xref: /petsc/src/mat/impls/sell/seq/fdsell.c (revision 3ba1676111f5c958fe6c2729b46ca4d523958bb3)
1d4002b98SHong Zhang #include <../src/mat/impls/sell/seq/sell.h>
2d4002b98SHong Zhang #include <../src/mat/impls/aij/seq/aij.h>
3d4002b98SHong Zhang #include <petsc/private/isimpl.h>
4d4002b98SHong Zhang 
5d4002b98SHong Zhang /*
6d4002b98SHong Zhang  MatGetColumnIJ_SeqSELL_Color() and MatRestoreColumnIJ_SeqSELL_Color() are customized from
7d4002b98SHong Zhang  MatGetColumnIJ_SeqSELL() and MatRestoreColumnIJ_SeqSELL() by adding an output
8d4002b98SHong Zhang  spidx[], index of a->a, to be used in MatTransposeColoringCreate_SeqSELL() and MatFDColoringCreate_SeqSELL()
9d4002b98SHong Zhang */
10d71ae5a4SJacob Faibussowitsch PetscErrorCode MatGetColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *nn, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done)
11d71ae5a4SJacob Faibussowitsch {
12d4002b98SHong Zhang   Mat_SeqSELL *a = (Mat_SeqSELL *)A->data;
13d4002b98SHong Zhang   PetscInt     i, j, *collengths, *cia, *cja, n = A->cmap->n, totalslices;
14d4002b98SHong Zhang   PetscInt     row, col;
15d4002b98SHong Zhang   PetscInt    *cspidx;
16d4002b98SHong Zhang   PetscBool    isnonzero;
17d4002b98SHong Zhang 
18d4002b98SHong Zhang   PetscFunctionBegin;
19d4002b98SHong Zhang   *nn = n;
20*3ba16761SJacob Faibussowitsch   if (!ia) PetscFunctionReturn(PETSC_SUCCESS);
21d4002b98SHong Zhang 
229566063dSJacob Faibussowitsch   PetscCall(PetscCalloc1(n + 1, &collengths));
239566063dSJacob Faibussowitsch   PetscCall(PetscMalloc1(n + 1, &cia));
249566063dSJacob Faibussowitsch   PetscCall(PetscMalloc1(a->nz + 1, &cja));
259566063dSJacob Faibussowitsch   PetscCall(PetscMalloc1(a->nz + 1, &cspidx));
26d4002b98SHong Zhang 
27d4002b98SHong Zhang   totalslices = A->rmap->n / 8 + ((A->rmap->n & 0x07) ? 1 : 0); /* floor(n/8) */
28d4002b98SHong Zhang   for (i = 0; i < totalslices; i++) {                           /* loop over slices */
29d4002b98SHong Zhang     for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
30d4002b98SHong Zhang       isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
31d4002b98SHong Zhang       if (isnonzero) collengths[a->colidx[j]]++;
32d4002b98SHong Zhang     }
33d4002b98SHong Zhang   }
34d4002b98SHong Zhang 
35d4002b98SHong Zhang   cia[0] = oshift;
36ad540459SPierre Jolivet   for (i = 0; i < n; i++) cia[i + 1] = cia[i] + collengths[i];
379566063dSJacob Faibussowitsch   PetscCall(PetscArrayzero(collengths, n));
38d4002b98SHong Zhang 
39d4002b98SHong Zhang   for (i = 0; i < totalslices; i++) { /* loop over slices */
40d4002b98SHong Zhang     for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) {
41d4002b98SHong Zhang       isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]);
42d4002b98SHong Zhang       if (isnonzero) {
43d4002b98SHong Zhang         col                                         = a->colidx[j];
44d4002b98SHong Zhang         cspidx[cia[col] + collengths[col] - oshift] = j;                    /* index of a->colidx */
45d4002b98SHong Zhang         cja[cia[col] + collengths[col] - oshift]    = 8 * i + row + oshift; /* row index */
46d4002b98SHong Zhang         collengths[col]++;
47d4002b98SHong Zhang       }
48d4002b98SHong Zhang     }
49d4002b98SHong Zhang   }
50d4002b98SHong Zhang 
519566063dSJacob Faibussowitsch   PetscCall(PetscFree(collengths));
529371c9d4SSatish Balay   *ia    = cia;
539371c9d4SSatish Balay   *ja    = cja;
54d4002b98SHong Zhang   *spidx = cspidx;
55*3ba16761SJacob Faibussowitsch   PetscFunctionReturn(PETSC_SUCCESS);
56d4002b98SHong Zhang }
57d4002b98SHong Zhang 
58d71ae5a4SJacob Faibussowitsch PetscErrorCode MatRestoreColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *n, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done)
59d71ae5a4SJacob Faibussowitsch {
60d4002b98SHong Zhang   PetscFunctionBegin;
61d4002b98SHong Zhang 
62*3ba16761SJacob Faibussowitsch   if (!ia) PetscFunctionReturn(PETSC_SUCCESS);
639566063dSJacob Faibussowitsch   PetscCall(PetscFree(*ia));
649566063dSJacob Faibussowitsch   PetscCall(PetscFree(*ja));
659566063dSJacob Faibussowitsch   PetscCall(PetscFree(*spidx));
66*3ba16761SJacob Faibussowitsch   PetscFunctionReturn(PETSC_SUCCESS);
67d4002b98SHong Zhang }
68