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_Inode(Mat,IS*,IS*);
7: EXTERN PetscErrorCode MatInodeGetInodeSizes_Inode(Mat,PetscInt*,PetscInt*[],PetscInt*);
12: PetscErrorCode MatView_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_Inode(Mat A, MatAssemblyType mode)
38: {
40: PetscTruth samestructure;
43: /* info.nz_unneeded of zero denotes no structural change was made to the matrix during Assembly */
44: samestructure = (PetscTruth)(!A->info.nz_unneeded);
45: /* check for identical nodes. If found, use inode functions */
46: Mat_CheckInode(A,samestructure);
47: return(0);
48: }
52: PetscErrorCode MatDestroy_Inode(Mat A)
53: {
55: Mat_SeqAIJ *a=(Mat_SeqAIJ*)A->data;
58: PetscFree(a->inode.size);
59: PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeAdjustForInodes_C","",PETSC_NULL);
60: PetscObjectComposeFunctionDynamic((PetscObject)A,"MatInodeGetInodeSizes_C","",PETSC_NULL);
61: return(0);
62: }
64: /* MatCreate_Inode is not DLLEXPORTed because it is not a constructor for a complete type. */
65: /* It is also not registered as a type for use within MatSetType. */
66: /* It is intended as a helper for the MATSEQAIJ class, so classes which desire Inodes should */
67: /* inherit off of MATSEQAIJ instead by calling MatSetType(MATSEQAIJ) in their constructor. */
68: /* Maybe this is a bad idea. (?) */
71: PetscErrorCode MatCreate_Inode(Mat B)
72: {
73: Mat_SeqAIJ *b=(Mat_SeqAIJ*)B->data;
77: b->inode.use = PETSC_TRUE;
78: b->inode.node_count = 0;
79: b->inode.size = 0;
80: b->inode.limit = 5;
81: b->inode.max_limit = 5;
83: PetscOptionsBegin(B->comm,B->prefix,"Options for SEQAIJ matrix","Mat");
84: PetscOptionsTruth("-mat_no_unroll","Do not optimize for inodes (slower)",PETSC_NULL,b->inode.use,&b->inode.use,PETSC_NULL);
85: if (!b->inode.use) {PetscInfo(B,"Not using Inode routines due to -mat_no_unroll\n");}
86: PetscOptionsTruth("-mat_no_inode","Do not optimize for inodes (slower)",PETSC_NULL,b->inode.use,&b->inode.use,PETSC_NULL);
87: if (!b->inode.use) {PetscInfo(B,"Not using Inode routines due to -mat_no_inode\n");}
88: PetscOptionsInt("-mat_inode_limit","Do not use inodes larger then this value",PETSC_NULL,b->inode.limit,&b->inode.limit,PETSC_NULL);
89: PetscOptionsEnd();
90: if (b->inode.limit > b->inode.max_limit) b->inode.limit = b->inode.max_limit;
92: PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeAdjustForInodes_C",
93: "MatInodeAdjustForInodes_Inode",
94: MatInodeAdjustForInodes_Inode);
95: PetscObjectComposeFunctionDynamic((PetscObject)B,"MatInodeGetInodeSizes_C",
96: "MatInodeGetInodeSizes_Inode",
97: MatInodeGetInodeSizes_Inode);
98: return(0);
99: }
103: PetscErrorCode MatSetOption_Inode(Mat A,MatOption op)
104: {
105: Mat_SeqAIJ *a=(Mat_SeqAIJ*)A->data;
109: switch(op) {
110: case MAT_USE_INODES:
111: a->inode.use = PETSC_TRUE;
112: break;
113: case MAT_DO_NOT_USE_INODES:
114: a->inode.use = PETSC_FALSE;
115: PetscInfo(A,"Not using Inode routines due to MatSetOption(MAT_DO_NOT_USE_INODES\n");
116: break;
117: case MAT_INODE_LIMIT_1:
118: a->inode.limit = 1;
119: break;
120: case MAT_INODE_LIMIT_2:
121: a->inode.limit = 2;
122: break;
123: case MAT_INODE_LIMIT_3:
124: a->inode.limit = 3;
125: break;
126: case MAT_INODE_LIMIT_4:
127: a->inode.limit = 4;
128: break;
129: case MAT_INODE_LIMIT_5:
130: a->inode.limit = 5;
131: break;
132: default:
133: break;
134: }
135: return(0);
136: }
140: PetscErrorCode MatDuplicate_Inode(Mat A,MatDuplicateOption cpvalues,Mat *C)
141: {
142: Mat B=*C;
143: Mat_SeqAIJ *c=(Mat_SeqAIJ*)B->data,*a=(Mat_SeqAIJ*)A->data;
145: PetscInt m=A->rmap.n;
149: c->inode.use = a->inode.use;
150: c->inode.limit = a->inode.limit;
151: c->inode.max_limit = a->inode.max_limit;
152: if (a->inode.size){
153: PetscMalloc((m+1)*sizeof(PetscInt),&c->inode.size);
154: c->inode.node_count = a->inode.node_count;
155: PetscMemcpy(c->inode.size,a->inode.size,(m+1)*sizeof(PetscInt));
156: } else {
157: c->inode.size = 0;
158: c->inode.node_count = 0;
159: }
160: return(0);
161: }
165: PetscErrorCode MatILUDTFactor_Inode(Mat A,IS isrow,IS iscol,MatFactorInfo *info,Mat *fact)
166: {
170: /* check for identical nodes. If found, use inode functions */
171: Mat_CheckInode(*fact,PETSC_FALSE);
172: return(0);
173: }
177: PetscErrorCode MatLUFactorSymbolic_Inode(Mat A,IS isrow,IS iscol,MatFactorInfo *info,Mat *fact)
178: {
182: /* check for identical nodes. If found, use inode functions */
183: Mat_CheckInode(*fact,PETSC_FALSE);
184: return(0);
185: }
189: PetscErrorCode MatILUFactorSymbolic_Inode(Mat A,IS isrow,IS iscol,MatFactorInfo *info,Mat *fact)
190: {
194: /* check for identical nodes. If found, use inode functions */
195: Mat_CheckInode(*fact,PETSC_FALSE);
196: return(0);
197: }