1 2 #include "private/matimpl.h" /*I "petscmat.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "MatCheckCompressedRow" 6 /*@C 7 MatCheckCompressedRow - Determines whether the compressed row matrix format should be used. 8 If the format is to be used, this routine creates Mat_CompressedRow struct. 9 Compressed row format provides high performance routines by taking advantage of zero rows. 10 Supported types are MATAIJ, MATBAIJ and MATSBAIJ. 11 12 Collective 13 14 Input Parameters: 15 + A - the matrix 16 . compressedrow - pointer to the struct Mat_CompressedRow 17 . ai - row pointer used by seqaij and seqbaij 18 . mbs - number of (block) rows represented by ai 19 - ratio - ratio of (num of zero rows)/m, used to determine if the compressed row format should be used 20 21 Notes: By default PETSc will not check for compressed rows on sequential matrices. Call MatSetOption(Mat,MAT_CHECK_COMPRESSED_ROW,PETSC_TRUE); before 22 MatAssemblyBegin() to have it check. 23 24 Level: developer 25 @*/ 26 PetscErrorCode MatCheckCompressedRow(Mat A,Mat_CompressedRow *compressedrow,PetscInt *ai,PetscInt mbs,PetscReal ratio) 27 { 28 PetscErrorCode ierr; 29 PetscInt nrows,*cpi=PETSC_NULL,*ridx=PETSC_NULL,nz,i,row; 30 31 PetscFunctionBegin; 32 if (!compressedrow->check) PetscFunctionReturn(0); 33 34 /* in case this is being reused, delete old space */ 35 ierr = PetscFree2(compressedrow->i,compressedrow->rindex);CHKERRQ(ierr); 36 compressedrow->i = PETSC_NULL; 37 compressedrow->rindex = PETSC_NULL; 38 39 40 /* compute number of zero rows */ 41 nrows = 0; 42 for (i=0; i<mbs; i++){ /* for each row */ 43 nz = ai[i+1] - ai[i]; /* number of nonzeros */ 44 if (nz == 0) nrows++; 45 } 46 /* if a large number of zero rows is found, use compressedrow data structure */ 47 if (nrows < ratio*mbs) { 48 compressedrow->use = PETSC_FALSE; 49 ierr = PetscInfo3(A,"Found the ratio (num_zerorows %d)/(num_localrows %d) < %G. Do not use CompressedRow routines.\n",nrows,mbs,ratio);CHKERRQ(ierr); 50 } else { 51 compressedrow->use = PETSC_TRUE; 52 ierr = PetscInfo3(A,"Found the ratio (num_zerorows %d)/(num_localrows %d) > %G. Use CompressedRow routines.\n",nrows,mbs,ratio);CHKERRQ(ierr); 53 54 /* set compressed row format */ 55 nrows = mbs - nrows; /* num of non-zero rows */ 56 ierr = PetscMalloc2(nrows+1,PetscInt,&cpi,nrows,PetscInt,&ridx);CHKERRQ(ierr); 57 row = 0; 58 cpi[0] = 0; 59 for (i=0; i<mbs; i++){ 60 nz = ai[i+1] - ai[i]; 61 if (nz == 0) continue; 62 cpi[row+1] = ai[i+1]; /* compressed row pointer */ 63 ridx[row++] = i; /* compressed row local index */ 64 } 65 compressedrow->nrows = nrows; 66 compressedrow->i = cpi; 67 compressedrow->rindex = ridx; 68 } 69 70 PetscFunctionReturn(0); 71 } 72