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