12a6744ebSBarry Smith /*
22a6744ebSBarry Smith This provides a matrix that applies a VecScatter to a vector.
32a6744ebSBarry Smith */
42a6744ebSBarry Smith
5af0996ceSBarry Smith #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/
6af0996ceSBarry Smith #include <petsc/private/vecimpl.h>
72a6744ebSBarry Smith
82a6744ebSBarry Smith typedef struct {
92a6744ebSBarry Smith VecScatter scatter;
102a6744ebSBarry Smith } Mat_Scatter;
112a6744ebSBarry Smith
122a6744ebSBarry Smith /*@
1311a5261eSBarry Smith MatScatterGetVecScatter - Returns the user-provided scatter set with `MatScatterSetVecScatter()` in a `MATSCATTER` matrix
142a6744ebSBarry Smith
152ef1f0ffSBarry Smith Logically Collective
162a6744ebSBarry Smith
172a6744ebSBarry Smith Input Parameter:
1811a5261eSBarry Smith . mat - the matrix, should have been created with MatCreateScatter() or have type `MATSCATTER`
192a6744ebSBarry Smith
202a6744ebSBarry Smith Output Parameter:
212a6744ebSBarry Smith . scatter - the scatter context
222a6744ebSBarry Smith
232a6744ebSBarry Smith Level: intermediate
242a6744ebSBarry Smith
25fe59aa6dSJacob Faibussowitsch .seealso: [](ch_matrices), `Mat`, `MATSCATTER`, `MatCreateScatter()`, `MatScatterSetVecScatter()`
262a6744ebSBarry Smith @*/
MatScatterGetVecScatter(Mat mat,VecScatter * scatter)27d71ae5a4SJacob Faibussowitsch PetscErrorCode MatScatterGetVecScatter(Mat mat, VecScatter *scatter)
28d71ae5a4SJacob Faibussowitsch {
292a6744ebSBarry Smith Mat_Scatter *mscatter;
302a6744ebSBarry Smith
312a6744ebSBarry Smith PetscFunctionBegin;
320700a824SBarry Smith PetscValidHeaderSpecific(mat, MAT_CLASSID, 1);
334f572ea9SToby Isaac PetscAssertPointer(scatter, 2);
342a6744ebSBarry Smith mscatter = (Mat_Scatter *)mat->data;
352a6744ebSBarry Smith *scatter = mscatter->scatter;
363ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
372a6744ebSBarry Smith }
382a6744ebSBarry Smith
MatDestroy_Scatter(Mat mat)3966976f2fSJacob Faibussowitsch static PetscErrorCode MatDestroy_Scatter(Mat mat)
40d71ae5a4SJacob Faibussowitsch {
412a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter *)mat->data;
422a6744ebSBarry Smith
432a6744ebSBarry Smith PetscFunctionBegin;
449566063dSJacob Faibussowitsch PetscCall(VecScatterDestroy(&scatter->scatter));
459566063dSJacob Faibussowitsch PetscCall(PetscFree(mat->data));
463ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
472a6744ebSBarry Smith }
482a6744ebSBarry Smith
MatMult_Scatter(Mat A,Vec x,Vec y)4966976f2fSJacob Faibussowitsch static PetscErrorCode MatMult_Scatter(Mat A, Vec x, Vec y)
50d71ae5a4SJacob Faibussowitsch {
512a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter *)A->data;
522a6744ebSBarry Smith
532a6744ebSBarry Smith PetscFunctionBegin;
5428b400f6SJacob Faibussowitsch PetscCheck(scatter->scatter, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Need to first call MatScatterSetScatter()");
559566063dSJacob Faibussowitsch PetscCall(VecZeroEntries(y));
569566063dSJacob Faibussowitsch PetscCall(VecScatterBegin(scatter->scatter, x, y, ADD_VALUES, SCATTER_FORWARD));
579566063dSJacob Faibussowitsch PetscCall(VecScatterEnd(scatter->scatter, x, y, ADD_VALUES, SCATTER_FORWARD));
583ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
592a6744ebSBarry Smith }
602a6744ebSBarry Smith
MatMultAdd_Scatter(Mat A,Vec x,Vec y,Vec z)6166976f2fSJacob Faibussowitsch static PetscErrorCode MatMultAdd_Scatter(Mat A, Vec x, Vec y, Vec z)
62d71ae5a4SJacob Faibussowitsch {
632a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter *)A->data;
642a6744ebSBarry Smith
652a6744ebSBarry Smith PetscFunctionBegin;
6628b400f6SJacob Faibussowitsch PetscCheck(scatter->scatter, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Need to first call MatScatterSetScatter()");
679566063dSJacob Faibussowitsch if (z != y) PetscCall(VecCopy(y, z));
689566063dSJacob Faibussowitsch PetscCall(VecScatterBegin(scatter->scatter, x, z, ADD_VALUES, SCATTER_FORWARD));
699566063dSJacob Faibussowitsch PetscCall(VecScatterEnd(scatter->scatter, x, z, ADD_VALUES, SCATTER_FORWARD));
703ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
712a6744ebSBarry Smith }
722a6744ebSBarry Smith
MatMultTranspose_Scatter(Mat A,Vec x,Vec y)7366976f2fSJacob Faibussowitsch static PetscErrorCode MatMultTranspose_Scatter(Mat A, Vec x, Vec y)
74d71ae5a4SJacob Faibussowitsch {
752a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter *)A->data;
762a6744ebSBarry Smith
772a6744ebSBarry Smith PetscFunctionBegin;
7828b400f6SJacob Faibussowitsch PetscCheck(scatter->scatter, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Need to first call MatScatterSetScatter()");
799566063dSJacob Faibussowitsch PetscCall(VecZeroEntries(y));
809566063dSJacob Faibussowitsch PetscCall(VecScatterBegin(scatter->scatter, x, y, ADD_VALUES, SCATTER_REVERSE));
819566063dSJacob Faibussowitsch PetscCall(VecScatterEnd(scatter->scatter, x, y, ADD_VALUES, SCATTER_REVERSE));
823ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
832a6744ebSBarry Smith }
842a6744ebSBarry Smith
MatMultTransposeAdd_Scatter(Mat A,Vec x,Vec y,Vec z)8566976f2fSJacob Faibussowitsch static PetscErrorCode MatMultTransposeAdd_Scatter(Mat A, Vec x, Vec y, Vec z)
86d71ae5a4SJacob Faibussowitsch {
872a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter *)A->data;
882a6744ebSBarry Smith
892a6744ebSBarry Smith PetscFunctionBegin;
9028b400f6SJacob Faibussowitsch PetscCheck(scatter->scatter, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "Need to first call MatScatterSetScatter()");
919566063dSJacob Faibussowitsch if (z != y) PetscCall(VecCopy(y, z));
929566063dSJacob Faibussowitsch PetscCall(VecScatterBegin(scatter->scatter, x, z, ADD_VALUES, SCATTER_REVERSE));
939566063dSJacob Faibussowitsch PetscCall(VecScatterEnd(scatter->scatter, x, z, ADD_VALUES, SCATTER_REVERSE));
943ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
952a6744ebSBarry Smith }
962a6744ebSBarry Smith
97f4259b30SLisandro Dalcin static struct _MatOps MatOps_Values = {NULL,
98f4259b30SLisandro Dalcin NULL,
99f4259b30SLisandro Dalcin NULL,
1002a6744ebSBarry Smith MatMult_Scatter,
1012a6744ebSBarry Smith /* 4*/ MatMultAdd_Scatter,
1022a6744ebSBarry Smith MatMultTranspose_Scatter,
1032a6744ebSBarry Smith MatMultTransposeAdd_Scatter,
104f4259b30SLisandro Dalcin NULL,
105f4259b30SLisandro Dalcin NULL,
106f4259b30SLisandro Dalcin NULL,
107f4259b30SLisandro Dalcin /* 10*/ NULL,
108f4259b30SLisandro Dalcin NULL,
109f4259b30SLisandro Dalcin NULL,
110f4259b30SLisandro Dalcin NULL,
111f4259b30SLisandro Dalcin NULL,
112f4259b30SLisandro Dalcin /* 15*/ NULL,
113f4259b30SLisandro Dalcin NULL,
114f4259b30SLisandro Dalcin NULL,
115f4259b30SLisandro Dalcin NULL,
116f4259b30SLisandro Dalcin NULL,
117f4259b30SLisandro Dalcin /* 20*/ NULL,
118f4259b30SLisandro Dalcin NULL,
119f4259b30SLisandro Dalcin NULL,
120f4259b30SLisandro Dalcin NULL,
121f4259b30SLisandro Dalcin /* 24*/ NULL,
122f4259b30SLisandro Dalcin NULL,
123f4259b30SLisandro Dalcin NULL,
124f4259b30SLisandro Dalcin NULL,
125f4259b30SLisandro Dalcin NULL,
126f4259b30SLisandro Dalcin /* 29*/ NULL,
127f4259b30SLisandro Dalcin NULL,
128f4259b30SLisandro Dalcin NULL,
129f4259b30SLisandro Dalcin NULL,
130f4259b30SLisandro Dalcin NULL,
131f4259b30SLisandro Dalcin /* 34*/ NULL,
132f4259b30SLisandro Dalcin NULL,
133f4259b30SLisandro Dalcin NULL,
134f4259b30SLisandro Dalcin NULL,
135f4259b30SLisandro Dalcin NULL,
136f4259b30SLisandro Dalcin /* 39*/ NULL,
137f4259b30SLisandro Dalcin NULL,
138f4259b30SLisandro Dalcin NULL,
139f4259b30SLisandro Dalcin NULL,
140f4259b30SLisandro Dalcin NULL,
141f4259b30SLisandro Dalcin /* 44*/ NULL,
142f4259b30SLisandro Dalcin NULL,
1437d68702bSBarry Smith MatShift_Basic,
144f4259b30SLisandro Dalcin NULL,
145f4259b30SLisandro Dalcin NULL,
146f4259b30SLisandro Dalcin /* 49*/ NULL,
147f4259b30SLisandro Dalcin NULL,
148f4259b30SLisandro Dalcin NULL,
149f4259b30SLisandro Dalcin NULL,
150f4259b30SLisandro Dalcin NULL,
151f4259b30SLisandro Dalcin /* 54*/ NULL,
152f4259b30SLisandro Dalcin NULL,
153f4259b30SLisandro Dalcin NULL,
154f4259b30SLisandro Dalcin NULL,
155f4259b30SLisandro Dalcin NULL,
156f4259b30SLisandro Dalcin /* 59*/ NULL,
1572a6744ebSBarry Smith MatDestroy_Scatter,
158f4259b30SLisandro Dalcin NULL,
159f4259b30SLisandro Dalcin NULL,
160f4259b30SLisandro Dalcin NULL,
161f4259b30SLisandro Dalcin /* 64*/ NULL,
162f4259b30SLisandro Dalcin NULL,
163f4259b30SLisandro Dalcin NULL,
164f4259b30SLisandro Dalcin NULL,
165f4259b30SLisandro Dalcin NULL,
166f4259b30SLisandro Dalcin /* 69*/ NULL,
167f4259b30SLisandro Dalcin NULL,
168f4259b30SLisandro Dalcin NULL,
169f4259b30SLisandro Dalcin NULL,
170f4259b30SLisandro Dalcin NULL,
171f4259b30SLisandro Dalcin /* 74*/ NULL,
172f4259b30SLisandro Dalcin NULL,
173f4259b30SLisandro Dalcin NULL,
174f4259b30SLisandro Dalcin NULL,
175f4259b30SLisandro Dalcin NULL,
176f4259b30SLisandro Dalcin /* 79*/ NULL,
177f4259b30SLisandro Dalcin NULL,
178f4259b30SLisandro Dalcin NULL,
179f4259b30SLisandro Dalcin NULL,
180f4259b30SLisandro Dalcin NULL,
181f4259b30SLisandro Dalcin /* 84*/ NULL,
182f4259b30SLisandro Dalcin NULL,
183f4259b30SLisandro Dalcin NULL,
184f4259b30SLisandro Dalcin NULL,
185f4259b30SLisandro Dalcin NULL,
186f4259b30SLisandro Dalcin /* 89*/ NULL,
187f4259b30SLisandro Dalcin NULL,
188f4259b30SLisandro Dalcin NULL,
189f4259b30SLisandro Dalcin NULL,
190f4259b30SLisandro Dalcin NULL,
191f4259b30SLisandro Dalcin /* 94*/ NULL,
192f4259b30SLisandro Dalcin NULL,
193f4259b30SLisandro Dalcin NULL,
194f4259b30SLisandro Dalcin NULL,
195f4259b30SLisandro Dalcin NULL,
196f4259b30SLisandro Dalcin /*99*/ NULL,
197f4259b30SLisandro Dalcin NULL,
198f4259b30SLisandro Dalcin NULL,
199f4259b30SLisandro Dalcin NULL,
200f4259b30SLisandro Dalcin NULL,
201f4259b30SLisandro Dalcin /*104*/ NULL,
202f4259b30SLisandro Dalcin NULL,
203f4259b30SLisandro Dalcin NULL,
204f4259b30SLisandro Dalcin NULL,
205f4259b30SLisandro Dalcin NULL,
206f4259b30SLisandro Dalcin /*109*/ NULL,
207f4259b30SLisandro Dalcin NULL,
208f4259b30SLisandro Dalcin NULL,
209f4259b30SLisandro Dalcin NULL,
210f4259b30SLisandro Dalcin NULL,
211f4259b30SLisandro Dalcin /*114*/ NULL,
212f4259b30SLisandro Dalcin NULL,
213f4259b30SLisandro Dalcin NULL,
214f4259b30SLisandro Dalcin NULL,
215f4259b30SLisandro Dalcin NULL,
216f4259b30SLisandro Dalcin /*119*/ NULL,
217f4259b30SLisandro Dalcin NULL,
218f4259b30SLisandro Dalcin NULL,
219f4259b30SLisandro Dalcin NULL,
220f4259b30SLisandro Dalcin NULL,
221f4259b30SLisandro Dalcin /*124*/ NULL,
222f4259b30SLisandro Dalcin NULL,
223f4259b30SLisandro Dalcin NULL,
224f4259b30SLisandro Dalcin NULL,
225f4259b30SLisandro Dalcin NULL,
226f4259b30SLisandro Dalcin /*129*/ NULL,
227f4259b30SLisandro Dalcin NULL,
228f4259b30SLisandro Dalcin NULL,
229f4259b30SLisandro Dalcin NULL,
230f4259b30SLisandro Dalcin NULL,
231f4259b30SLisandro Dalcin /*134*/ NULL,
232f4259b30SLisandro Dalcin NULL,
233f4259b30SLisandro Dalcin NULL,
234f4259b30SLisandro Dalcin NULL,
235f4259b30SLisandro Dalcin NULL,
236f4259b30SLisandro Dalcin /*139*/ NULL,
237d70f29a3SPierre Jolivet NULL,
238d70f29a3SPierre Jolivet NULL,
239*03db1824SAlex Lindsay NULL,
240dec0b466SHong Zhang NULL};
2412a6744ebSBarry Smith
2422a6744ebSBarry Smith /*MC
2432ef1f0ffSBarry Smith MATSCATTER - "scatter" - A matrix type that simply applies a `VecScatterBegin()` and `VecScatterEnd()` to perform `MatMult()`
2442a6744ebSBarry Smith
2452a6744ebSBarry Smith Level: advanced
2462a6744ebSBarry Smith
247283d8c04SPierre Jolivet .seealso: [](ch_matrices), `Mat`, `MATSCATTER`, MatCreateScatter()`, `MatScatterSetVecScatter()`, `MatScatterGetVecScatter()`
2482a6744ebSBarry Smith M*/
2492a6744ebSBarry Smith
MatCreate_Scatter(Mat A)250d71ae5a4SJacob Faibussowitsch PETSC_EXTERN PetscErrorCode MatCreate_Scatter(Mat A)
251d71ae5a4SJacob Faibussowitsch {
2522a6744ebSBarry Smith Mat_Scatter *b;
2532a6744ebSBarry Smith
2542a6744ebSBarry Smith PetscFunctionBegin;
255aea10558SJacob Faibussowitsch A->ops[0] = MatOps_Values;
2564dfa11a4SJacob Faibussowitsch PetscCall(PetscNew(&b));
2572a6744ebSBarry Smith
2582a6744ebSBarry Smith A->data = (void *)b;
2592a6744ebSBarry Smith
2609566063dSJacob Faibussowitsch PetscCall(PetscLayoutSetUp(A->rmap));
2619566063dSJacob Faibussowitsch PetscCall(PetscLayoutSetUp(A->cmap));
2622a6744ebSBarry Smith
2632a6744ebSBarry Smith A->assembled = PETSC_TRUE;
2642a6744ebSBarry Smith A->preallocated = PETSC_FALSE;
26574d95942SJed Brown
2669566063dSJacob Faibussowitsch PetscCall(PetscObjectChangeTypeName((PetscObject)A, MATSCATTER));
2673ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
2682a6744ebSBarry Smith }
2692a6744ebSBarry Smith
27097929ea7SJunchao Zhang #include <petsc/private/sfimpl.h>
271cc4c1da9SBarry Smith /*@
27211a5261eSBarry Smith MatCreateScatter - Creates a new matrix of `MatType` `MATSCATTER`, based on a VecScatter
2732a6744ebSBarry Smith
274d083f849SBarry Smith Collective
2752a6744ebSBarry Smith
2762a6744ebSBarry Smith Input Parameters:
2772a6744ebSBarry Smith + comm - MPI communicator
27811a5261eSBarry Smith - scatter - a `VecScatter`
2792a6744ebSBarry Smith
2802a6744ebSBarry Smith Output Parameter:
2812a6744ebSBarry Smith . A - the matrix
2822a6744ebSBarry Smith
2832a6744ebSBarry Smith Level: intermediate
2842a6744ebSBarry Smith
285a3b724e8SBarry Smith Notes:
2862a6744ebSBarry Smith PETSc requires that matrices and vectors being used for certain
2872a6744ebSBarry Smith operations are partitioned accordingly. For example, when
2882a6744ebSBarry Smith creating a scatter matrix, A, that supports parallel matrix-vector
28911a5261eSBarry Smith products using `MatMult`(A,x,y) the user should set the number
2902a6744ebSBarry Smith of local matrix rows to be the number of local elements of the
2912a6744ebSBarry Smith corresponding result vector, y. Note that this is information is
2922a6744ebSBarry Smith required for use of the matrix interface routines, even though
2932a6744ebSBarry Smith the scatter matrix may not actually be physically partitioned.
2942a6744ebSBarry Smith
295fe59aa6dSJacob Faibussowitsch Developer Notes:
29611a5261eSBarry Smith This directly accesses information inside the `VecScatter` associated with the matrix-vector product
2976eb45d04SBarry Smith for this matrix. This is not desirable..
2986eb45d04SBarry Smith
2991cc06b55SBarry Smith .seealso: [](ch_matrices), `Mat`, `MatScatterSetVecScatter()`, `MatScatterGetVecScatter()`, `MATSCATTER`
3002a6744ebSBarry Smith @*/
MatCreateScatter(MPI_Comm comm,VecScatter scatter,Mat * A)301d71ae5a4SJacob Faibussowitsch PetscErrorCode MatCreateScatter(MPI_Comm comm, VecScatter scatter, Mat *A)
302d71ae5a4SJacob Faibussowitsch {
3032a6744ebSBarry Smith PetscFunctionBegin;
3049566063dSJacob Faibussowitsch PetscCall(MatCreate(comm, A));
3059566063dSJacob Faibussowitsch PetscCall(MatSetSizes(*A, scatter->vscat.to_n, scatter->vscat.from_n, PETSC_DETERMINE, PETSC_DETERMINE));
3069566063dSJacob Faibussowitsch PetscCall(MatSetType(*A, MATSCATTER));
3079566063dSJacob Faibussowitsch PetscCall(MatScatterSetVecScatter(*A, scatter));
3089566063dSJacob Faibussowitsch PetscCall(MatSetUp(*A));
3093ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
3102a6744ebSBarry Smith }
3112a6744ebSBarry Smith
3122a6744ebSBarry Smith /*@
31311a5261eSBarry Smith MatScatterSetVecScatter - sets the scatter that the matrix is to apply as its linear operator in a `MATSCATTER`
3142a6744ebSBarry Smith
3152ef1f0ffSBarry Smith Logically Collective
3162a6744ebSBarry Smith
3172a6744ebSBarry Smith Input Parameters:
31811a5261eSBarry Smith + mat - the `MATSCATTER` matrix
31911a5261eSBarry Smith - scatter - the scatter context create with `VecScatterCreate()`
3202a6744ebSBarry Smith
3212a6744ebSBarry Smith Level: advanced
3222a6744ebSBarry Smith
323fe59aa6dSJacob Faibussowitsch .seealso: [](ch_matrices), `Mat`, `MATSCATTER`, `MatCreateScatter()`
3242a6744ebSBarry Smith @*/
MatScatterSetVecScatter(Mat mat,VecScatter scatter)325d71ae5a4SJacob Faibussowitsch PetscErrorCode MatScatterSetVecScatter(Mat mat, VecScatter scatter)
326d71ae5a4SJacob Faibussowitsch {
3272a6744ebSBarry Smith Mat_Scatter *mscatter = (Mat_Scatter *)mat->data;
3282a6744ebSBarry Smith
3292a6744ebSBarry Smith PetscFunctionBegin;
3300700a824SBarry Smith PetscValidHeaderSpecific(mat, MAT_CLASSID, 1);
33197929ea7SJunchao Zhang PetscValidHeaderSpecific(scatter, PETSCSF_CLASSID, 2);
332835f2295SStefano Zampini PetscCheckSameComm(scatter, 2, mat, 1);
33308401ef6SPierre Jolivet PetscCheck(mat->rmap->n == scatter->vscat.to_n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Number of local rows in matrix %" PetscInt_FMT " not equal local scatter size %" PetscInt_FMT, mat->rmap->n, scatter->vscat.to_n);
33408401ef6SPierre Jolivet PetscCheck(mat->cmap->n == scatter->vscat.from_n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Number of local columns in matrix %" PetscInt_FMT " not equal local scatter size %" PetscInt_FMT, mat->cmap->n, scatter->vscat.from_n);
3352a6744ebSBarry Smith
3369566063dSJacob Faibussowitsch PetscCall(PetscObjectReference((PetscObject)scatter));
3379566063dSJacob Faibussowitsch PetscCall(VecScatterDestroy(&mscatter->scatter));
33826fbe8dcSKarl Rupp
3392a6744ebSBarry Smith mscatter->scatter = scatter;
3403ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
3412a6744ebSBarry Smith }
342