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