xref: /petsc/src/mat/impls/aij/mpi/aijsell/mpiaijsell.c (revision 9c5460f9064ca60dd71a234a1f6faf93e7a6b0c9)
1 #include <../src/mat/impls/aij/mpi/mpiaij.h>
2 /*@C
3   MatCreateMPIAIJSELL - Creates a sparse parallel matrix whose local
4   portions are stored as `MATSEQAIJSELL` matrices (a matrix class that inherits
5   from SEQAIJ but performs some operations in SELL format).
6 
7   Collective
8 
9   Input Parameters:
10 + comm  - MPI communicator
11 . m     - number of local rows (or `PETSC_DECIDE` to have calculated if `M` is given)
12            This value should be the same as the local size used in creating the
13            y vector for the matrix-vector product y = Ax.
14 . n     - This value should be the same as the local size used in creating the
15        x vector for the matrix-vector product y = Ax. (or `PETSC_DECIDE` to have
16        calculated if `N` is given) For square matrices `n` is almost always `m`.
17 . M     - number of global rows (or `PETSC_DETERMINE` to have calculated if `m` is given)
18 . N     - number of global columns (or `PETSC_DETERMINE` to have calculated if `n` is given)
19 . d_nz  - number of nonzeros per row in DIAGONAL portion of local submatrix
20            (same value is used for all local rows)
21 . d_nnz - array containing the number of nonzeros in the various rows of the
22            DIAGONAL portion of the local submatrix (possibly different for each row)
23            or `NULL`, if `d_nz` is used to specify the nonzero structure.
24            The size of this array is equal to the number of local rows, i.e `m`.
25            For matrices you plan to factor you must leave room for the diagonal entry and
26            put in the entry even if it is zero.
27 . o_nz  - number of nonzeros per row in the OFF-DIAGONAL portion of local
28            submatrix (same value is used for all local rows).
29 - o_nnz - array containing the number of nonzeros in the various rows of the
30            OFF-DIAGONAL portion of the local submatrix (possibly different for
31            each row) or `NULL`, if `o_nz` is used to specify the nonzero
32            structure. The size of this array is equal to the number
33            of local rows, i.e `m`.
34 
35   Output Parameter:
36 . A - the matrix
37 
38   Options Database Key:
39 . -mat_aijsell_eager_shadow - Construct shadow matrix upon matrix assembly; default is to take a "lazy" approach, performing this step the first
40                                time the matrix is applied
41 
42   Level: intermediate
43 
44   Notes:
45   If the *_nnz parameter is given then the *_nz parameter is ignored
46 
47   `m`,`n`,`M`,`N` parameters specify the size of the matrix, and its partitioning across
48   processors, while `d_nz`,`d_nnz`,`o_nz`,`o_nnz` parameters specify the approximate
49   storage requirements for this matrix.
50 
51   If `PETSC_DECIDE` or `PETSC_DETERMINE` is used for a particular argument on one
52   processor than it must be used on all processors that share the object for
53   that argument.
54 
55   The user MUST specify either the local or global matrix dimensions
56   (possibly both).
57 
58   The parallel matrix is partitioned such that the first m0 rows belong to
59   process 0, the next m1 rows belong to process 1, the next m2 rows belong
60   to process 2 etc.. where m0,m1,m2... are the input parameter `m`.
61 
62   The DIAGONAL portion of the local submatrix of a processor can be defined
63   as the submatrix which is obtained by extraction the part corresponding
64   to the rows r1-r2 and columns r1-r2 of the global matrix, where r1 is the
65   first row that belongs to the processor, and r2 is the last row belonging
66   to the this processor. This is a square mxm matrix. The remaining portion
67   of the local submatrix (mxN) constitute the OFF-DIAGONAL portion.
68 
69   If `o_nnz`, `d_nnz` are specified, then `o_nz`, and `d_nz` are ignored.
70 
71   When calling this routine with a single process communicator, a matrix of
72   type `MATSEQAIJSELL` is returned.  If a matrix of type `MATMPIAIJSELL` is desired
73   for this type of communicator, use the construction mechanism
74 .vb
75    MatCreate(...,&A);
76    MatSetType(A,MPIAIJSELL);
77    MatMPIAIJSetPreallocation(A,...);
78 .ve
79 
80 .seealso: [](ch_matrices), `Mat`, [Sparse Matrix Creation](sec_matsparse), `MATSEQAIJSELL`, `MATMPIAIJSELL`, `MATAIJSELL`, `MatCreate()`, `MatCreateSeqAIJSELL()`, `MatSetValues()`
81 @*/
MatCreateMPIAIJSELL(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[],Mat * A)82 PetscErrorCode MatCreateMPIAIJSELL(MPI_Comm comm, PetscInt m, PetscInt n, PetscInt M, PetscInt N, PetscInt d_nz, const PetscInt d_nnz[], PetscInt o_nz, const PetscInt o_nnz[], Mat *A)
83 {
84   PetscMPIInt size;
85 
86   PetscFunctionBegin;
87   PetscCall(MatCreate(comm, A));
88   PetscCall(MatSetSizes(*A, m, n, M, N));
89   PetscCallMPI(MPI_Comm_size(comm, &size));
90   if (size > 1) {
91     PetscCall(MatSetType(*A, MATMPIAIJSELL));
92     PetscCall(MatMPIAIJSetPreallocation(*A, d_nz, d_nnz, o_nz, o_nnz));
93   } else {
94     PetscCall(MatSetType(*A, MATSEQAIJSELL));
95     PetscCall(MatSeqAIJSetPreallocation(*A, d_nz, d_nnz));
96   }
97   PetscFunctionReturn(PETSC_SUCCESS);
98 }
99 
100 PETSC_INTERN PetscErrorCode MatConvert_SeqAIJ_SeqAIJSELL(Mat, MatType, MatReuse, Mat *);
101 
MatMPIAIJSetPreallocation_MPIAIJSELL(Mat B,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[])102 static PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJSELL(Mat B, PetscInt d_nz, const PetscInt d_nnz[], PetscInt o_nz, const PetscInt o_nnz[])
103 {
104   Mat_MPIAIJ *b = (Mat_MPIAIJ *)B->data;
105 
106   PetscFunctionBegin;
107   PetscCall(MatMPIAIJSetPreallocation_MPIAIJ(B, d_nz, d_nnz, o_nz, o_nnz));
108   PetscCall(MatConvert_SeqAIJ_SeqAIJSELL(b->A, MATSEQAIJSELL, MAT_INPLACE_MATRIX, &b->A));
109   PetscCall(MatConvert_SeqAIJ_SeqAIJSELL(b->B, MATSEQAIJSELL, MAT_INPLACE_MATRIX, &b->B));
110   PetscFunctionReturn(PETSC_SUCCESS);
111 }
112 
MatConvert_MPIAIJ_MPIAIJSELL(Mat A,MatType type,MatReuse reuse,Mat * newmat)113 PETSC_INTERN PetscErrorCode MatConvert_MPIAIJ_MPIAIJSELL(Mat A, MatType type, MatReuse reuse, Mat *newmat)
114 {
115   Mat B = *newmat;
116 
117   PetscFunctionBegin;
118   if (reuse == MAT_INITIAL_MATRIX) PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &B));
119 
120   PetscCall(PetscObjectChangeTypeName((PetscObject)B, MATMPIAIJSELL));
121   PetscCall(PetscObjectComposeFunction((PetscObject)B, "MatMPIAIJSetPreallocation_C", MatMPIAIJSetPreallocation_MPIAIJSELL));
122   *newmat = B;
123   PetscFunctionReturn(PETSC_SUCCESS);
124 }
125 
MatCreate_MPIAIJSELL(Mat A)126 PETSC_EXTERN PetscErrorCode MatCreate_MPIAIJSELL(Mat A)
127 {
128   PetscFunctionBegin;
129   PetscCall(MatSetType(A, MATMPIAIJ));
130   PetscCall(MatConvert_MPIAIJ_MPIAIJSELL(A, MATMPIAIJSELL, MAT_INPLACE_MATRIX, &A));
131   PetscFunctionReturn(PETSC_SUCCESS);
132 }
133 
134 /*MC
135    MATAIJSELL - "AIJSELL" - A matrix type to be used for sparse matrices.
136 
137    This matrix type is identical to `MATSEQAIJSELL` when constructed with a single process communicator,
138    and `MATMPIAIJSELL` otherwise.  As a result, for single process communicators,
139    MatSeqAIJSetPreallocation() is supported, and similarly `MatMPIAIJSetPreallocation()` is supported
140    for communicators controlling multiple processes.  It is recommended that you call both of
141    the above preallocation routines for simplicity.
142 
143    Options Database Key:
144 . -mat_type aijsell - sets the matrix type to `MATAIJSELL`
145 
146   Level: beginner
147 
148 .seealso: [](ch_matrices), `Mat`, `MatCreateMPIAIJSELL()`, `MATSEQAIJSELL`, `MATMPIAIJSELL`, `MATSEQAIJ`, `MATMPIAIJ`, `MATSEQAIJPERM`, `MATMPIAIJPERM`, `MATSEQAIJMKL`, `MATMPIAIJMKL`
149 M*/
150