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 *nn = n; 20 if (!ia) PetscFunctionReturn(0); 21 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)); 26 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 } 34 35 cia[0] = oshift; 36 for (i=0; i<n; i++) { 37 cia[i+1] = cia[i] + collengths[i]; 38 } 39 PetscCall(PetscArrayzero(collengths,n)); 40 41 for (i=0; i<totalslices; i++) { /* loop over slices */ 42 for (j=a->sliidx[i],row=0; j<a->sliidx[i+1]; j++,row=((row+1)&0x07)) { 43 isnonzero = (PetscBool)((j-a->sliidx[i])/8 < a->rlen[8*i+row]); 44 if (isnonzero) { 45 col = a->colidx[j]; 46 cspidx[cia[col]+collengths[col]-oshift] = j; /* index of a->colidx */ 47 cja[cia[col]+collengths[col]-oshift] = 8*i+row +oshift; /* row index */ 48 collengths[col]++; 49 } 50 } 51 } 52 53 PetscCall(PetscFree(collengths)); 54 *ia = cia; *ja = cja; 55 *spidx = cspidx; 56 PetscFunctionReturn(0); 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 63 if (!ia) PetscFunctionReturn(0); 64 PetscCall(PetscFree(*ia)); 65 PetscCall(PetscFree(*ja)); 66 PetscCall(PetscFree(*spidx)); 67 PetscFunctionReturn(0); 68 } 69