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