Actual source code: inode2.c

  1: #define PETSCMAT_DLL
 2:  #include ../src/mat/impls/aij/seq/aij.h

  4: EXTERN PetscErrorCode Mat_CheckInode(Mat,PetscTruth);
  6: EXTERN PetscErrorCode  MatInodeAdjustForInodes_SeqAIJ_Inode(Mat,IS*,IS*);
  7: EXTERN PetscErrorCode  MatInodeGetInodeSizes_SeqAIJ_Inode(Mat,PetscInt*,PetscInt*[],PetscInt*);

 12: PetscErrorCode MatView_SeqAIJ_Inode(Mat A,PetscViewer viewer)
 13: {
 14:   Mat_SeqAIJ         *a=(Mat_SeqAIJ*)A->data;
 15:   PetscErrorCode    ierr;
 16:   PetscTruth        iascii;
 17:   PetscViewerFormat format;

 20:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
 21:   if (iascii) {
 22:     PetscViewerGetFormat(viewer,&format);
 23:     if (format == PETSC_VIEWER_ASCII_INFO_DETAIL || format == PETSC_VIEWER_ASCII_INFO) {
 24:       if (a->inode.size) {
 25:         PetscViewerASCIIPrintf(viewer,"using I-node routines: found %D nodes, limit used is %D\n",
 26:                                       a->inode.node_count,a->inode.limit);
 27:       } else {
 28:         PetscViewerASCIIPrintf(viewer,"not using I-node routines\n");
 29:       }
 30:     }
 31:   }
 32:   return(0);
 33: }

 37: PetscErrorCode MatAssemblyEnd_SeqAIJ_Inode(Mat A, MatAssemblyType mode)
 38: {
 39:   Mat_SeqAIJ     *a = (Mat_SeqAIJ*)A->data;
 41:   PetscTruth     samestructure;

 44:   /* info.nz_unneeded of zero denotes no structural change was made to the matrix during Assembly */
 45:   samestructure = (PetscTruth)(!A->info.nz_unneeded);
 46:   /* check for identical nodes. If found, use inode functions */
 47:   Mat_CheckInode(A,samestructure);
 48:   a->inode.ibdiagvalid = PETSC_FALSE;
 49:   return(0);
 50: }

 54: PetscErrorCode MatDestroy_SeqAIJ_Inode(Mat A)
 55: {
 57:   Mat_SeqAIJ     *a=(Mat_SeqAIJ*)A->data;

 60:   PetscFree(a->inode.size);
 61:   PetscFree3(a->inode.ibdiag,a->inode.bdiag,a->inode.ssor_work);
 62:   PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeAdjustForInodes_C","",PETSC_NULL);
 63:   PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeGetInodeSizes_C","",PETSC_NULL);
 64:   return(0);
 65: }

 67: /* MatCreate_SeqAIJ_Inode is not DLLEXPORTed because it is not a constructor for a complete type.    */
 68: /* It is also not registered as a type for use within MatSetType.                             */
 69: /* It is intended as a helper for the MATSEQAIJ class, so classes which desire Inodes should  */
 70: /*    inherit off of MATSEQAIJ instead by calling MatSetType(MATSEQAIJ) in their constructor. */
 71: /* Maybe this is a bad idea. (?) */
 74: PetscErrorCode MatCreate_SeqAIJ_Inode(Mat B)
 75: {
 76:   Mat_SeqAIJ     *b=(Mat_SeqAIJ*)B->data;
 78:   PetscTruth     no_inode,no_unroll;

 81:   no_inode             = PETSC_FALSE;
 82:   no_unroll            = PETSC_FALSE;
 83:   b->inode.node_count  = 0;
 84:   b->inode.size        = 0;
 85:   b->inode.limit       = 5;
 86:   b->inode.max_limit   = 5;
 87:   b->inode.ibdiagvalid = PETSC_FALSE;
 88:   b->inode.ibdiag      = 0;
 89:   b->inode.bdiag       = 0;

 91:   PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Options for SEQAIJ matrix","Mat");
 92:     PetscOptionsTruth("-mat_no_unroll","Do not optimize for inodes (slower)",PETSC_NULL,no_unroll,&no_unroll,PETSC_NULL);
 93:     if (no_unroll) {PetscInfo(B,"Not using Inode routines due to -mat_no_unroll\n");}
 94:     PetscOptionsTruth("-mat_no_inode","Do not optimize for inodes (slower)",PETSC_NULL,no_inode,&no_inode,PETSC_NULL);
 95:     if (no_inode) {PetscInfo(B,"Not using Inode routines due to -mat_no_inode\n");}
 96:     PetscOptionsInt("-mat_inode_limit","Do not use inodes larger then this value",PETSC_NULL,b->inode.limit,&b->inode.limit,PETSC_NULL);
 97:   PetscOptionsEnd();
 98:   b->inode.use = (PetscTruth)(!(no_unroll || no_inode));
 99:   if (b->inode.limit > b->inode.max_limit) b->inode.limit = b->inode.max_limit;

101:   PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeAdjustForInodes_C",
102:                                      "MatInodeAdjustForInodes_SeqAIJ_Inode",
103:                                       MatInodeAdjustForInodes_SeqAIJ_Inode);
104:   PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeGetInodeSizes_C",
105:                                      "MatInodeGetInodeSizes_SeqAIJ_Inode",
106:                                       MatInodeGetInodeSizes_SeqAIJ_Inode);
107:   return(0);
108: }

112: PetscErrorCode MatSetOption_SeqAIJ_Inode(Mat A,MatOption op,PetscTruth flg)
113: {
114:   Mat_SeqAIJ     *a=(Mat_SeqAIJ*)A->data;

117:   switch(op) {
118:     case MAT_USE_INODES:
119:       a->inode.use         = flg;
120:       break;
121:     default:
122:       break;
123:   }
124:   return(0);
125: }

129: PetscErrorCode MatDuplicate_SeqAIJ_Inode(Mat A,MatDuplicateOption cpvalues,Mat *C)
130: {
131:   Mat            B=*C;
132:   Mat_SeqAIJ     *c=(Mat_SeqAIJ*)B->data,*a=(Mat_SeqAIJ*)A->data;
134:   PetscInt       m=A->rmap->n;


138:   c->inode.use          = a->inode.use;
139:   c->inode.limit        = a->inode.limit;
140:   c->inode.max_limit    = a->inode.max_limit;
141:   if (a->inode.size){
142:     PetscMalloc((m+1)*sizeof(PetscInt),&c->inode.size);
143:     c->inode.node_count = a->inode.node_count;
144:     PetscMemcpy(c->inode.size,a->inode.size,(m+1)*sizeof(PetscInt));
145:   } else {
146:     c->inode.size       = 0;
147:     c->inode.node_count = 0;
148:   }
149:   c->inode.ibdiagvalid = PETSC_FALSE;
150:   c->inode.ibdiag      = 0;
151:   c->inode.bdiag       = 0;
152:   return(0);
153: }