1 #include <../src/mat/impls/sell/seq/sell.h> 2 #include <../src/mat/impls/aij/seq/aij.h> 3 #include <petsc/private/isimpl.h> 4 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; 17 18 PetscFunctionBegin; 19 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Expected cmap->n %" PetscInt_FMT " >= 0", n); 20 *nn = n; 21 if (!ia) PetscFunctionReturn(PETSC_SUCCESS); 22 23 PetscCall(PetscCalloc1(n, &collengths)); 24 PetscCall(PetscMalloc1(n + 1, &cia)); 25 PetscCall(PetscMalloc1(a->nz + 1, &cja)); 26 PetscCall(PetscMalloc1(a->nz + 1, &cspidx)); 27 28 totalslices = PetscCeilInt(A->rmap->n, 8); 29 for (i = 0; i < totalslices; i++) { /* loop over slices */ 30 for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) { 31 isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]); 32 if (isnonzero) collengths[a->colidx[j]]++; 33 } 34 } 35 36 cia[0] = oshift; 37 for (i = 0; i < n; i++) cia[i + 1] = cia[i] + collengths[i]; 38 PetscCall(PetscArrayzero(collengths, n)); 39 40 for (i = 0; i < totalslices; i++) { /* loop over slices */ 41 for (j = a->sliidx[i], row = 0; j < a->sliidx[i + 1]; j++, row = ((row + 1) & 0x07)) { 42 isnonzero = (PetscBool)((j - a->sliidx[i]) / 8 < a->rlen[8 * i + row]); 43 if (isnonzero) { 44 col = a->colidx[j]; 45 cspidx[cia[col] + collengths[col] - oshift] = j; /* index of a->colidx */ 46 cja[cia[col] + collengths[col] - oshift] = 8 * i + row + oshift; /* row index */ 47 collengths[col]++; 48 } 49 } 50 } 51 52 PetscCall(PetscFree(collengths)); 53 *ia = cia; 54 *ja = cja; 55 *spidx = cspidx; 56 PetscFunctionReturn(PETSC_SUCCESS); 57 } 58 59 PetscErrorCode MatRestoreColumnIJ_SeqSELL_Color(Mat A, PetscInt oshift, PetscBool symmetric, PetscBool inodecompressed, PetscInt *n, const PetscInt *ia[], const PetscInt *ja[], PetscInt *spidx[], PetscBool *done) 60 { 61 PetscFunctionBegin; 62 if (!ia) PetscFunctionReturn(PETSC_SUCCESS); 63 PetscCall(PetscFree(*ia)); 64 PetscCall(PetscFree(*ja)); 65 PetscCall(PetscFree(*spidx)); 66 PetscFunctionReturn(PETSC_SUCCESS); 67 } 68