1*d4002b98SHong Zhang #include <../src/mat/impls/sell/seq/sell.h> 2*d4002b98SHong Zhang #include <../src/mat/impls/aij/seq/aij.h> 3*d4002b98SHong Zhang #include <petsc/private/isimpl.h> 4*d4002b98SHong Zhang 5*d4002b98SHong Zhang /* 6*d4002b98SHong Zhang MatGetColumnIJ_SeqSELL_Color() and MatRestoreColumnIJ_SeqSELL_Color() are customized from 7*d4002b98SHong Zhang MatGetColumnIJ_SeqSELL() and MatRestoreColumnIJ_SeqSELL() by adding an output 8*d4002b98SHong Zhang spidx[], index of a->a, to be used in MatTransposeColoringCreate_SeqSELL() and MatFDColoringCreate_SeqSELL() 9*d4002b98SHong Zhang */ 10*d4002b98SHong Zhang 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*d4002b98SHong Zhang { 12*d4002b98SHong Zhang Mat_SeqSELL *a = (Mat_SeqSELL*)A->data; 13*d4002b98SHong Zhang PetscInt i,j,*collengths,*cia,*cja,n = A->cmap->n,totalslices; 14*d4002b98SHong Zhang PetscInt row,col; 15*d4002b98SHong Zhang PetscInt *cspidx; 16*d4002b98SHong Zhang PetscBool isnonzero; 17*d4002b98SHong Zhang PetscErrorCode ierr; 18*d4002b98SHong Zhang 19*d4002b98SHong Zhang PetscFunctionBegin; 20*d4002b98SHong Zhang *nn = n; 21*d4002b98SHong Zhang if (!ia) PetscFunctionReturn(0); 22*d4002b98SHong Zhang 23*d4002b98SHong Zhang ierr = PetscCalloc1(n+1,&collengths);CHKERRQ(ierr); 24*d4002b98SHong Zhang ierr = PetscMalloc1(n+1,&cia);CHKERRQ(ierr); 25*d4002b98SHong Zhang ierr = PetscMalloc1(a->nz+1,&cja);CHKERRQ(ierr); 26*d4002b98SHong Zhang ierr = PetscMalloc1(a->nz+1,&cspidx);CHKERRQ(ierr); 27*d4002b98SHong Zhang 28*d4002b98SHong Zhang totalslices = A->rmap->n/8+((A->rmap->n & 0x07)?1:0); /* floor(n/8) */ 29*d4002b98SHong Zhang for (i=0; i<totalslices; i++) { /* loop over slices */ 30*d4002b98SHong Zhang for (j=a->sliidx[i],row=0; j<a->sliidx[i+1]; j++,row=((row+1)&0x07)) { 31*d4002b98SHong Zhang isnonzero = (PetscBool)((j-a->sliidx[i])/8 < a->rlen[8*i+row]); 32*d4002b98SHong Zhang if (isnonzero) collengths[a->colidx[j]]++; 33*d4002b98SHong Zhang } 34*d4002b98SHong Zhang } 35*d4002b98SHong Zhang 36*d4002b98SHong Zhang cia[0] = oshift; 37*d4002b98SHong Zhang for (i=0; i<n; i++) { 38*d4002b98SHong Zhang cia[i+1] = cia[i] + collengths[i]; 39*d4002b98SHong Zhang } 40*d4002b98SHong Zhang ierr = PetscMemzero(collengths,n*sizeof(PetscInt));CHKERRQ(ierr); 41*d4002b98SHong Zhang 42*d4002b98SHong Zhang for (i=0; i<totalslices; i++) { /* loop over slices */ 43*d4002b98SHong Zhang for (j=a->sliidx[i],row=0; j<a->sliidx[i+1]; j++,row=((row+1)&0x07)) { 44*d4002b98SHong Zhang isnonzero = (PetscBool)((j-a->sliidx[i])/8 < a->rlen[8*i+row]); 45*d4002b98SHong Zhang if (isnonzero) { 46*d4002b98SHong Zhang col = a->colidx[j]; 47*d4002b98SHong Zhang cspidx[cia[col]+collengths[col]-oshift] = j; /* index of a->colidx */ 48*d4002b98SHong Zhang cja[cia[col]+collengths[col]-oshift] = 8*i+row +oshift; /* row index */ 49*d4002b98SHong Zhang collengths[col]++; 50*d4002b98SHong Zhang } 51*d4002b98SHong Zhang } 52*d4002b98SHong Zhang } 53*d4002b98SHong Zhang 54*d4002b98SHong Zhang ierr = PetscFree(collengths);CHKERRQ(ierr); 55*d4002b98SHong Zhang *ia = cia; *ja = cja; 56*d4002b98SHong Zhang *spidx = cspidx; 57*d4002b98SHong Zhang PetscFunctionReturn(0); 58*d4002b98SHong Zhang } 59*d4002b98SHong Zhang 60*d4002b98SHong Zhang PetscErrorCode MatRestoreColumnIJ_SeqSELL_Color(Mat A,PetscInt oshift,PetscBool symmetric,PetscBool inodecompressed,PetscInt *n,const PetscInt *ia[],const PetscInt *ja[],PetscInt *spidx[],PetscBool *done) 61*d4002b98SHong Zhang { 62*d4002b98SHong Zhang PetscErrorCode ierr; 63*d4002b98SHong Zhang 64*d4002b98SHong Zhang PetscFunctionBegin; 65*d4002b98SHong Zhang 66*d4002b98SHong Zhang if (!ia) PetscFunctionReturn(0); 67*d4002b98SHong Zhang ierr = PetscFree(*ia);CHKERRQ(ierr); 68*d4002b98SHong Zhang ierr = PetscFree(*ja);CHKERRQ(ierr); 69*d4002b98SHong Zhang ierr = PetscFree(*spidx);CHKERRQ(ierr); 70*d4002b98SHong Zhang PetscFunctionReturn(0); 71*d4002b98SHong Zhang } 72