1 2 #include <../src/mat/impls/aij/seq/aij.h> 3 4 extern PetscErrorCode Mat_CheckInode(Mat,PetscBool); 5 EXTERN_C_BEGIN 6 extern PetscErrorCode MatInodeAdjustForInodes_SeqAIJ_Inode(Mat,IS*,IS*); 7 extern PetscErrorCode MatInodeGetInodeSizes_SeqAIJ_Inode(Mat,PetscInt*,PetscInt*[],PetscInt*); 8 EXTERN_C_END 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "MatView_SeqAIJ_Inode" 12 PetscErrorCode MatView_SeqAIJ_Inode(Mat A,PetscViewer viewer) 13 { 14 Mat_SeqAIJ *a=(Mat_SeqAIJ*)A->data; 15 PetscErrorCode ierr; 16 PetscBool iascii; 17 PetscViewerFormat format; 18 19 PetscFunctionBegin; 20 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 21 if (iascii) { 22 ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr); 23 if (format == PETSC_VIEWER_ASCII_INFO_DETAIL || format == PETSC_VIEWER_ASCII_INFO) { 24 if (a->inode.size) { 25 ierr = PetscViewerASCIIPrintf(viewer,"using I-node routines: found %D nodes, limit used is %D\n", 26 a->inode.node_count,a->inode.limit);CHKERRQ(ierr); 27 } else { 28 ierr = PetscViewerASCIIPrintf(viewer,"not using I-node routines\n");CHKERRQ(ierr); 29 } 30 } 31 } 32 PetscFunctionReturn(0); 33 } 34 35 #undef __FUNCT__ 36 #define __FUNCT__ "MatAssemblyEnd_SeqAIJ_Inode" 37 PetscErrorCode MatAssemblyEnd_SeqAIJ_Inode(Mat A, MatAssemblyType mode) 38 { 39 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 40 PetscErrorCode ierr; 41 PetscBool samestructure; 42 43 PetscFunctionBegin; 44 /* info.nz_unneeded of zero denotes no structural change was made to the matrix during Assembly */ 45 samestructure = (PetscBool)(!A->info.nz_unneeded); 46 /* check for identical nodes. If found, use inode functions */ 47 ierr = Mat_CheckInode(A,samestructure);CHKERRQ(ierr); 48 49 a->inode.ibdiagvalid = PETSC_FALSE; 50 PetscFunctionReturn(0); 51 } 52 53 #undef __FUNCT__ 54 #define __FUNCT__ "MatDestroy_SeqAIJ_Inode" 55 PetscErrorCode MatDestroy_SeqAIJ_Inode(Mat A) 56 { 57 PetscErrorCode ierr; 58 Mat_SeqAIJ *a=(Mat_SeqAIJ*)A->data; 59 60 PetscFunctionBegin; 61 ierr = PetscFree(a->inode.size);CHKERRQ(ierr); 62 ierr = PetscFree3(a->inode.ibdiag,a->inode.bdiag,a->inode.ssor_work);CHKERRQ(ierr); 63 ierr = PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeAdjustForInodes_C","",PETSC_NULL);CHKERRQ(ierr); 64 ierr = PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeGetInodeSizes_C","",PETSC_NULL);CHKERRQ(ierr); 65 PetscFunctionReturn(0); 66 } 67 68 /* MatCreate_SeqAIJ_Inode is not DLLEXPORTed because it is not a constructor for a complete type. */ 69 /* It is also not registered as a type for use within MatSetType. */ 70 /* It is intended as a helper for the MATSEQAIJ class, so classes which desire Inodes should */ 71 /* inherit off of MATSEQAIJ instead by calling MatSetType(MATSEQAIJ) in their constructor. */ 72 /* Maybe this is a bad idea. (?) */ 73 #undef __FUNCT__ 74 #define __FUNCT__ "MatCreate_SeqAIJ_Inode" 75 PetscErrorCode MatCreate_SeqAIJ_Inode(Mat B) 76 { 77 Mat_SeqAIJ *b=(Mat_SeqAIJ*)B->data; 78 PetscErrorCode ierr; 79 PetscBool no_inode,no_unroll; 80 81 PetscFunctionBegin; 82 no_inode = PETSC_FALSE; 83 no_unroll = PETSC_FALSE; 84 b->inode.node_count = 0; 85 b->inode.size = 0; 86 b->inode.limit = 5; 87 b->inode.max_limit = 5; 88 b->inode.ibdiagvalid = PETSC_FALSE; 89 b->inode.ibdiag = 0; 90 b->inode.bdiag = 0; 91 92 ierr = PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Options for SEQAIJ matrix","Mat");CHKERRQ(ierr); 93 ierr = PetscOptionsBool("-mat_no_unroll","Do not optimize for inodes (slower)",PETSC_NULL,no_unroll,&no_unroll,PETSC_NULL);CHKERRQ(ierr); 94 if (no_unroll) { 95 ierr = PetscInfo(B,"Not using Inode routines due to -mat_no_unroll\n");CHKERRQ(ierr); 96 } 97 ierr = PetscOptionsBool("-mat_no_inode","Do not optimize for inodes -slower-",PETSC_NULL,no_inode,&no_inode,PETSC_NULL);CHKERRQ(ierr); 98 if (no_inode) { 99 ierr = PetscInfo(B,"Not using Inode routines due to -mat_no_inode\n");CHKERRQ(ierr); 100 } 101 ierr = PetscOptionsInt("-mat_inode_limit","Do not use inodes larger then this value",PETSC_NULL,b->inode.limit,&b->inode.limit,PETSC_NULL);CHKERRQ(ierr); 102 ierr = PetscOptionsEnd();CHKERRQ(ierr); 103 b->inode.use = (PetscBool)(!(no_unroll || no_inode)); 104 if (b->inode.limit > b->inode.max_limit) b->inode.limit = b->inode.max_limit; 105 106 ierr = PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeAdjustForInodes_C", 107 "MatInodeAdjustForInodes_SeqAIJ_Inode", 108 MatInodeAdjustForInodes_SeqAIJ_Inode);CHKERRQ(ierr); 109 ierr = PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeGetInodeSizes_C", 110 "MatInodeGetInodeSizes_SeqAIJ_Inode", 111 MatInodeGetInodeSizes_SeqAIJ_Inode);CHKERRQ(ierr); 112 PetscFunctionReturn(0); 113 } 114 115 #undef __FUNCT__ 116 #define __FUNCT__ "MatSetOption_SeqAIJ_Inode" 117 PetscErrorCode MatSetOption_SeqAIJ_Inode(Mat A,MatOption op,PetscBool flg) 118 { 119 Mat_SeqAIJ *a=(Mat_SeqAIJ*)A->data; 120 121 PetscFunctionBegin; 122 switch (op) { 123 case MAT_USE_INODES: 124 a->inode.use = flg; 125 break; 126 default: 127 break; 128 } 129 PetscFunctionReturn(0); 130 } 131 132 133 134 135 136