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: *nn = n;
19: if (!ia) return 0;
21: PetscCalloc1(n+1,&collengths);
22: PetscMalloc1(n+1,&cia);
23: PetscMalloc1(a->nz+1,&cja);
24: PetscMalloc1(a->nz+1,&cspidx);
26: totalslices = A->rmap->n/8+((A->rmap->n & 0x07)?1:0); /* floor(n/8) */
27: for (i=0; i<totalslices; i++) { /* loop over slices */
28: for (j=a->sliidx[i],row=0; j<a->sliidx[i+1]; j++,row=((row+1)&0x07)) {
29: isnonzero = (PetscBool)((j-a->sliidx[i])/8 < a->rlen[8*i+row]);
30: if (isnonzero) collengths[a->colidx[j]]++;
31: }
32: }
34: cia[0] = oshift;
35: for (i=0; i<n; i++) {
36: cia[i+1] = cia[i] + collengths[i];
37: }
38: PetscArrayzero(collengths,n);
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: }
52: PetscFree(collengths);
53: *ia = cia; *ja = cja;
54: *spidx = cspidx;
55: return 0;
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: {
61: if (!ia) return 0;
62: PetscFree(*ia);
63: PetscFree(*ja);
64: PetscFree(*spidx);
65: return 0;
66: }