xref: /petsc/src/mat/impls/composite/mcomposite.c (revision 04d534cea0490f06ff1ca11a5e66368ab98e561f)
1 
2 #include <petsc/private/matimpl.h>        /*I "petscmat.h" I*/
3 
4 typedef struct _Mat_CompositeLink *Mat_CompositeLink;
5 struct _Mat_CompositeLink {
6   Mat               mat;
7   Vec               work;
8   Mat_CompositeLink next,prev;
9 };
10 
11 typedef struct {
12   MatCompositeType  type;
13   Mat_CompositeLink head,tail;
14   Vec               work;
15   PetscScalar       scale;        /* scale factor supplied with MatScale() */
16   Vec               left,right;   /* left and right diagonal scaling provided with MatDiagonalScale() */
17   Vec               leftwork,rightwork;
18   PetscInt          nmat;
19   PetscBool         mergefromright;
20 } Mat_Composite;
21 
22 PetscErrorCode MatDestroy_Composite(Mat mat)
23 {
24   PetscErrorCode    ierr;
25   Mat_Composite     *shell = (Mat_Composite*)mat->data;
26   Mat_CompositeLink next   = shell->head,oldnext;
27 
28   PetscFunctionBegin;
29   while (next) {
30     ierr = MatDestroy(&next->mat);CHKERRQ(ierr);
31     if (next->work && (!next->next || next->work != next->next->work)) {
32       ierr = VecDestroy(&next->work);CHKERRQ(ierr);
33     }
34     oldnext = next;
35     next    = next->next;
36     ierr    = PetscFree(oldnext);CHKERRQ(ierr);
37   }
38   ierr = VecDestroy(&shell->work);CHKERRQ(ierr);
39   ierr = VecDestroy(&shell->left);CHKERRQ(ierr);
40   ierr = VecDestroy(&shell->right);CHKERRQ(ierr);
41   ierr = VecDestroy(&shell->leftwork);CHKERRQ(ierr);
42   ierr = VecDestroy(&shell->rightwork);CHKERRQ(ierr);
43   ierr = PetscFree(mat->data);CHKERRQ(ierr);
44   PetscFunctionReturn(0);
45 }
46 
47 PetscErrorCode MatMult_Composite_Multiplicative(Mat A,Vec x,Vec y)
48 {
49   Mat_Composite     *shell = (Mat_Composite*)A->data;
50   Mat_CompositeLink next   = shell->head;
51   PetscErrorCode    ierr;
52   Vec               in,out;
53 
54   PetscFunctionBegin;
55   if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
56   in = x;
57   if (shell->right) {
58     if (!shell->rightwork) {
59       ierr = VecDuplicate(shell->right,&shell->rightwork);CHKERRQ(ierr);
60     }
61     ierr = VecPointwiseMult(shell->rightwork,shell->right,in);CHKERRQ(ierr);
62     in   = shell->rightwork;
63   }
64   while (next->next) {
65     if (!next->work) { /* should reuse previous work if the same size */
66       ierr = MatCreateVecs(next->mat,NULL,&next->work);CHKERRQ(ierr);
67     }
68     out  = next->work;
69     ierr = MatMult(next->mat,in,out);CHKERRQ(ierr);
70     in   = out;
71     next = next->next;
72   }
73   ierr = MatMult(next->mat,in,y);CHKERRQ(ierr);
74   if (shell->left) {
75     ierr = VecPointwiseMult(y,shell->left,y);CHKERRQ(ierr);
76   }
77   ierr = VecScale(y,shell->scale);CHKERRQ(ierr);
78   PetscFunctionReturn(0);
79 }
80 
81 PetscErrorCode MatMultTranspose_Composite_Multiplicative(Mat A,Vec x,Vec y)
82 {
83   Mat_Composite     *shell = (Mat_Composite*)A->data;
84   Mat_CompositeLink tail   = shell->tail;
85   PetscErrorCode    ierr;
86   Vec               in,out;
87 
88   PetscFunctionBegin;
89   if (!tail) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
90   in = x;
91   if (shell->left) {
92     if (!shell->leftwork) {
93       ierr = VecDuplicate(shell->left,&shell->leftwork);CHKERRQ(ierr);
94     }
95     ierr = VecPointwiseMult(shell->leftwork,shell->left,in);CHKERRQ(ierr);
96     in   = shell->leftwork;
97   }
98   while (tail->prev) {
99     if (!tail->prev->work) { /* should reuse previous work if the same size */
100       ierr = MatCreateVecs(tail->mat,NULL,&tail->prev->work);CHKERRQ(ierr);
101     }
102     out  = tail->prev->work;
103     ierr = MatMultTranspose(tail->mat,in,out);CHKERRQ(ierr);
104     in   = out;
105     tail = tail->prev;
106   }
107   ierr = MatMultTranspose(tail->mat,in,y);CHKERRQ(ierr);
108   if (shell->right) {
109     ierr = VecPointwiseMult(y,shell->right,y);CHKERRQ(ierr);
110   }
111   ierr = VecScale(y,shell->scale);CHKERRQ(ierr);
112   PetscFunctionReturn(0);
113 }
114 
115 PetscErrorCode MatMult_Composite(Mat A,Vec x,Vec y)
116 {
117   Mat_Composite     *shell = (Mat_Composite*)A->data;
118   Mat_CompositeLink next   = shell->head;
119   PetscErrorCode    ierr;
120   Vec               in;
121 
122   PetscFunctionBegin;
123   if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
124   in = x;
125   if (shell->right) {
126     if (!shell->rightwork) {
127       ierr = VecDuplicate(shell->right,&shell->rightwork);CHKERRQ(ierr);
128     }
129     ierr = VecPointwiseMult(shell->rightwork,shell->right,in);CHKERRQ(ierr);
130     in   = shell->rightwork;
131   }
132   ierr = MatMult(next->mat,in,y);CHKERRQ(ierr);
133   while ((next = next->next)) {
134     ierr = MatMultAdd(next->mat,in,y,y);CHKERRQ(ierr);
135   }
136   if (shell->left) {
137     ierr = VecPointwiseMult(y,shell->left,y);CHKERRQ(ierr);
138   }
139   ierr = VecScale(y,shell->scale);CHKERRQ(ierr);
140   PetscFunctionReturn(0);
141 }
142 
143 PetscErrorCode MatMultTranspose_Composite(Mat A,Vec x,Vec y)
144 {
145   Mat_Composite     *shell = (Mat_Composite*)A->data;
146   Mat_CompositeLink next   = shell->head;
147   PetscErrorCode    ierr;
148   Vec               in;
149 
150   PetscFunctionBegin;
151   if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
152   in = x;
153   if (shell->left) {
154     if (!shell->leftwork) {
155       ierr = VecDuplicate(shell->left,&shell->leftwork);CHKERRQ(ierr);
156     }
157     ierr = VecPointwiseMult(shell->leftwork,shell->left,in);CHKERRQ(ierr);
158     in   = shell->leftwork;
159   }
160   ierr = MatMultTranspose(next->mat,in,y);CHKERRQ(ierr);
161   while ((next = next->next)) {
162     ierr = MatMultTransposeAdd(next->mat,in,y,y);CHKERRQ(ierr);
163   }
164   if (shell->right) {
165     ierr = VecPointwiseMult(y,shell->right,y);CHKERRQ(ierr);
166   }
167   ierr = VecScale(y,shell->scale);CHKERRQ(ierr);
168   PetscFunctionReturn(0);
169 }
170 
171 PetscErrorCode MatGetDiagonal_Composite(Mat A,Vec v)
172 {
173   Mat_Composite     *shell = (Mat_Composite*)A->data;
174   Mat_CompositeLink next   = shell->head;
175   PetscErrorCode    ierr;
176 
177   PetscFunctionBegin;
178   if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
179   if (shell->right || shell->left) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot get diagonal if left or right scaling");
180 
181   ierr = MatGetDiagonal(next->mat,v);CHKERRQ(ierr);
182   if (next->next && !shell->work) {
183     ierr = VecDuplicate(v,&shell->work);CHKERRQ(ierr);
184   }
185   while ((next = next->next)) {
186     ierr = MatGetDiagonal(next->mat,shell->work);CHKERRQ(ierr);
187     ierr = VecAXPY(v,1.0,shell->work);CHKERRQ(ierr);
188   }
189   ierr = VecScale(v,shell->scale);CHKERRQ(ierr);
190   PetscFunctionReturn(0);
191 }
192 
193 PetscErrorCode MatAssemblyEnd_Composite(Mat Y,MatAssemblyType t)
194 {
195   PetscErrorCode ierr;
196   PetscBool      flg = PETSC_FALSE;
197 
198   PetscFunctionBegin;
199   ierr = PetscOptionsGetBool(((PetscObject)Y)->options,((PetscObject)Y)->prefix,"-mat_composite_merge",&flg,NULL);CHKERRQ(ierr);
200   if (flg) {
201     ierr = MatCompositeMerge(Y);CHKERRQ(ierr);
202   }
203   PetscFunctionReturn(0);
204 }
205 
206 PetscErrorCode MatScale_Composite(Mat inA,PetscScalar alpha)
207 {
208   Mat_Composite *a = (Mat_Composite*)inA->data;
209 
210   PetscFunctionBegin;
211   a->scale *= alpha;
212   PetscFunctionReturn(0);
213 }
214 
215 PetscErrorCode MatDiagonalScale_Composite(Mat inA,Vec left,Vec right)
216 {
217   Mat_Composite  *a = (Mat_Composite*)inA->data;
218   PetscErrorCode ierr;
219 
220   PetscFunctionBegin;
221   if (left) {
222     if (!a->left) {
223       ierr = VecDuplicate(left,&a->left);CHKERRQ(ierr);
224       ierr = VecCopy(left,a->left);CHKERRQ(ierr);
225     } else {
226       ierr = VecPointwiseMult(a->left,left,a->left);CHKERRQ(ierr);
227     }
228   }
229   if (right) {
230     if (!a->right) {
231       ierr = VecDuplicate(right,&a->right);CHKERRQ(ierr);
232       ierr = VecCopy(right,a->right);CHKERRQ(ierr);
233     } else {
234       ierr = VecPointwiseMult(a->right,right,a->right);CHKERRQ(ierr);
235     }
236   }
237   PetscFunctionReturn(0);
238 }
239 
240 /*@
241    MatCreateComposite - Creates a matrix as the sum of zero or more matrices
242 
243   Collective on MPI_Comm
244 
245    Input Parameters:
246 +  comm - MPI communicator
247 .  nmat - number of matrices to put in
248 -  mats - the matrices
249 
250    Output Parameter:
251 .  mat - the matrix
252 
253    Level: advanced
254 
255    Notes:
256      Alternative construction
257 $       MatCreate(comm,&mat);
258 $       MatSetSizes(mat,m,n,M,N);
259 $       MatSetType(mat,MATCOMPOSITE);
260 $       MatCompositeAddMat(mat,mats[0]);
261 $       ....
262 $       MatCompositeAddMat(mat,mats[nmat-1]);
263 $       MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
264 $       MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
265 
266      For the multiplicative form the product is mat[nmat-1]*mat[nmat-2]*....*mat[0]
267 
268 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCompositeMerge(), MatCompositeSetType(), MatCompositeType
269 
270 @*/
271 PetscErrorCode MatCreateComposite(MPI_Comm comm,PetscInt nmat,const Mat *mats,Mat *mat)
272 {
273   PetscErrorCode ierr;
274   PetscInt       m,n,M,N,i;
275 
276   PetscFunctionBegin;
277   if (nmat < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must pass in at least one matrix");
278   PetscValidPointer(mat,3);
279 
280   ierr = MatGetLocalSize(mats[0],PETSC_IGNORE,&n);CHKERRQ(ierr);
281   ierr = MatGetLocalSize(mats[nmat-1],&m,PETSC_IGNORE);CHKERRQ(ierr);
282   ierr = MatGetSize(mats[0],PETSC_IGNORE,&N);CHKERRQ(ierr);
283   ierr = MatGetSize(mats[nmat-1],&M,PETSC_IGNORE);CHKERRQ(ierr);
284   ierr = MatCreate(comm,mat);CHKERRQ(ierr);
285   ierr = MatSetSizes(*mat,m,n,M,N);CHKERRQ(ierr);
286   ierr = MatSetType(*mat,MATCOMPOSITE);CHKERRQ(ierr);
287   for (i=0; i<nmat; i++) {
288     ierr = MatCompositeAddMat(*mat,mats[i]);CHKERRQ(ierr);
289   }
290   ierr = MatAssemblyBegin(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
291   ierr = MatAssemblyEnd(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
292   PetscFunctionReturn(0);
293 }
294 
295 
296 static PetscErrorCode MatCompositeAddMat_Composite(Mat mat,Mat smat)
297 {
298   Mat_Composite     *shell = (Mat_Composite*)mat->data;
299   Mat_CompositeLink ilink,next = shell->head;
300   PetscErrorCode    ierr;
301 
302   PetscFunctionBegin;
303   ierr        = PetscNewLog(mat,&ilink);CHKERRQ(ierr);
304   ilink->next = 0;
305   ierr        = PetscObjectReference((PetscObject)smat);CHKERRQ(ierr);
306   ilink->mat  = smat;
307 
308   if (!next) shell->head = ilink;
309   else {
310     while (next->next) {
311       next = next->next;
312     }
313     next->next  = ilink;
314     ilink->prev = next;
315   }
316   shell->tail =  ilink;
317   shell->nmat += 1;
318   PetscFunctionReturn(0);
319 }
320 
321 /*@
322     MatCompositeAddMat - add another matrix to a composite matrix
323 
324    Collective on Mat
325 
326     Input Parameters:
327 +   mat - the composite matrix
328 -   smat - the partial matrix
329 
330    Level: advanced
331 
332 .seealso: MatCreateComposite()
333 @*/
334 PetscErrorCode MatCompositeAddMat(Mat mat,Mat smat)
335 {
336   PetscErrorCode    ierr;
337 
338   PetscFunctionBegin;
339   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
340   PetscValidHeaderSpecific(smat,MAT_CLASSID,2);
341   ierr = PetscUseMethod(mat,"MatCompositeAddMat_C",(Mat,Mat),(mat,smat));CHKERRQ(ierr);
342   PetscFunctionReturn(0);
343 }
344 
345 static PetscErrorCode MatCompositeSetType_Composite(Mat mat,MatCompositeType type)
346 {
347   Mat_Composite  *b = (Mat_Composite*)mat->data;
348 
349   PetscFunctionBegin;
350   if (type == MAT_COMPOSITE_MULTIPLICATIVE) {
351     mat->ops->getdiagonal   = 0;
352     mat->ops->mult          = MatMult_Composite_Multiplicative;
353     mat->ops->multtranspose = MatMultTranspose_Composite_Multiplicative;
354     b->type                 = MAT_COMPOSITE_MULTIPLICATIVE;
355   } else {
356     mat->ops->getdiagonal   = MatGetDiagonal_Composite;
357     mat->ops->mult          = MatMult_Composite;
358     mat->ops->multtranspose = MatMultTranspose_Composite;
359     b->type                 = MAT_COMPOSITE_ADDITIVE;
360   }
361   PetscFunctionReturn(0);
362 }
363 
364 /*@
365    MatCompositeSetType - Indicates if the matrix is defined as the sum of a set of matrices or the product
366 
367   Collective on MPI_Comm
368 
369    Input Parameters:
370 .  mat - the composite matrix
371 
372 
373    Level: advanced
374 
375    Notes:
376       The MatType of the resulting matrix will be the same as the MatType of the FIRST
377     matrix in the composite matrix.
378 
379 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCreateComposite(), MATCOMPOSITE
380 
381 @*/
382 PetscErrorCode MatCompositeSetType(Mat mat,MatCompositeType type)
383 {
384   PetscErrorCode ierr;
385 
386   PetscFunctionBegin;
387   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
388   ierr = PetscUseMethod(mat,"MatCompositeSetType_C",(Mat,MatCompositeType),(mat,type));CHKERRQ(ierr);
389   PetscFunctionReturn(0);
390 }
391 
392 static PetscErrorCode MatCompositeMerge_Composite(Mat mat)
393 {
394   Mat_Composite     *shell = (Mat_Composite*)mat->data;
395   Mat_CompositeLink next   = shell->head, prev = shell->tail;
396   PetscErrorCode    ierr;
397   Mat               tmat,newmat;
398   Vec               left,right;
399   PetscScalar       scale;
400 
401   PetscFunctionBegin;
402   if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()");
403 
404   PetscFunctionBegin;
405   if (shell->type == MAT_COMPOSITE_ADDITIVE) {
406     ierr = MatDuplicate(next->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr);
407     while ((next = next->next)) {
408       ierr = MatAXPY(tmat,1.0,next->mat,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
409     }
410   } else {
411     if (shell->mergefromright) {
412       ierr = MatDuplicate(next->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr);
413       while ((next = next->next)) {
414         ierr = MatMatMult(next->mat,tmat,MAT_INITIAL_MATRIX,PETSC_DECIDE,&newmat);CHKERRQ(ierr);
415         ierr = MatDestroy(&tmat);CHKERRQ(ierr);
416         tmat = newmat;
417       }
418     } else {
419       ierr = MatDuplicate(prev->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr);
420       while ((prev = prev->prev)) {
421         ierr = MatMatMult(tmat,prev->mat,MAT_INITIAL_MATRIX,PETSC_DECIDE,&newmat);CHKERRQ(ierr);
422         ierr = MatDestroy(&tmat);CHKERRQ(ierr);
423         tmat = newmat;
424       }
425     }
426   }
427 
428   scale = shell->scale;
429   if ((left = shell->left)) {ierr = PetscObjectReference((PetscObject)left);CHKERRQ(ierr);}
430   if ((right = shell->right)) {ierr = PetscObjectReference((PetscObject)right);CHKERRQ(ierr);}
431 
432   ierr = MatHeaderReplace(mat,&tmat);CHKERRQ(ierr);
433 
434   ierr = MatDiagonalScale(mat,left,right);CHKERRQ(ierr);
435   ierr = MatScale(mat,scale);CHKERRQ(ierr);
436   ierr = VecDestroy(&left);CHKERRQ(ierr);
437   ierr = VecDestroy(&right);CHKERRQ(ierr);
438   PetscFunctionReturn(0);
439 }
440 
441 /*@
442    MatCompositeMerge - Given a composite matrix, replaces it with a "regular" matrix
443      by summing all the matrices inside the composite matrix.
444 
445   Collective on MPI_Comm
446 
447    Input Parameters:
448 .  mat - the composite matrix
449 
450 
451    Options Database:
452 .  -mat_composite_merge  (you must call MatAssemblyBegin()/MatAssemblyEnd() to have this checked)
453 
454    Level: advanced
455 
456    Notes:
457       The MatType of the resulting matrix will be the same as the MatType of the FIRST
458     matrix in the composite matrix.
459 
460 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCreateComposite(), MatCompositeSetMergeFromRight(), MATCOMPOSITE
461 
462 @*/
463 PetscErrorCode MatCompositeMerge(Mat mat)
464 {
465   PetscErrorCode ierr;
466 
467   PetscFunctionBegin;
468   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
469   ierr = PetscUseMethod(mat,"MatCompositeMerge_C",(Mat),(mat));CHKERRQ(ierr);
470   PetscFunctionReturn(0);
471 }
472 
473 static PetscErrorCode MatCompositeGetNMat_Composite(Mat mat,PetscInt *nmat)
474 {
475   Mat_Composite     *shell = (Mat_Composite*)mat->data;
476 
477   PetscFunctionBegin;
478   *nmat = shell->nmat;
479   PetscFunctionReturn(0);
480 }
481 
482 /*@
483    MatCompositeGetNMat - Returns the number of matrices in composite.
484 
485    Not Collective
486 
487    Input Parameter:
488 .  mat - the composite matrix
489 
490    Output Parameter:
491 .  size - the local size
492 .  nmat - number of matrices in composite
493 
494    Level: advanced
495 
496 .seealso: MatCreateComposite(), MatCompositeGetMat()
497 
498 @*/
499 PetscErrorCode MatCompositeGetNMat(Mat mat,PetscInt *nmat)
500 {
501   PetscErrorCode ierr;
502 
503   PetscFunctionBegin;
504   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
505   PetscValidPointer(nmat,2);
506   ierr = PetscUseMethod(mat,"MatCompositeGetNMat_C",(Mat,PetscInt*),(mat,nmat));CHKERRQ(ierr);
507   PetscFunctionReturn(0);
508 }
509 
510 static PetscErrorCode MatCompositeGetMat_Composite(Mat mat,PetscInt i,Mat *Ai)
511 {
512   Mat_Composite     *shell = (Mat_Composite*)mat->data;
513   Mat_CompositeLink ilink;
514   PetscInt          k;
515 
516   PetscFunctionBegin;
517   if (i >= shell->nmat) SETERRQ2(PetscObjectComm((PetscObject)mat),PETSC_ERR_ARG_OUTOFRANGE,"index out of range: %d >= %d",i,shell->nmat);
518   ilink = shell->head;
519   for (k=0; k<i; k++) {
520     ilink = ilink->next;
521   }
522   *Ai = ilink->mat;
523   PetscFunctionReturn(0);
524 }
525 
526 /*@
527    MatCompositeGetMat - Returns the ith matrix from composite.
528 
529    Logically Collective on Mat
530 
531    Input Parameter:
532 +  mat - the composite matrix
533 -  i - the number of requested matrix
534 
535    Output Parameter:
536 .  Ai - ith matrix in composite
537 
538    Level: advanced
539 
540 .seealso: MatCreateComposite(), MatCompositeGetNMat()
541 
542 @*/
543 PetscErrorCode MatCompositeGetMat(Mat mat,PetscInt i,Mat *Ai)
544 {
545   PetscErrorCode ierr;
546 
547   PetscFunctionBegin;
548   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
549   PetscValidLogicalCollectiveInt(mat,i,2);
550   PetscValidPointer(Ai,3);
551   ierr = PetscUseMethod(mat,"MatCompositeGetMat_C",(Mat,PetscInt,Mat*),(mat,i,Ai));CHKERRQ(ierr);
552   PetscFunctionReturn(0);
553 }
554 
555 static PetscErrorCode MatCompositeSetMergeFromRight_Composite(Mat mat,PetscBool flg)
556 {
557   Mat_Composite     *shell = (Mat_Composite*)mat->data;
558 
559   PetscFunctionBegin;
560   shell->mergefromright = flg;
561   PetscFunctionReturn(0);
562 }
563 
564 /*@
565    MatCompositeSetMergeFromRight - Sets MatCompositeMerge() to start from right.
566 
567    Logically Collective on Mat
568 
569    Input Parameter:
570 +  mat - the composite matrix
571 -  flg - if true (default) the merge starts with the first added matrix (mat[0])
572 
573    Level: advanced
574 
575 .seealso: MatCreateComposite(), MatCompositeMerge()
576 
577 @*/
578 PetscErrorCode MatCompositeSetMergeFromRight(Mat mat,PetscBool flg)
579 {
580   PetscErrorCode ierr;
581 
582   PetscFunctionBegin;
583   PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
584   PetscValidLogicalCollectiveBool(mat,flg,2);
585   ierr = PetscUseMethod(mat,"MatCompositeSetMergeFromRight_C",(Mat,PetscBool),(mat,flg));CHKERRQ(ierr);
586   PetscFunctionReturn(0);
587 }
588 
589 static struct _MatOps MatOps_Values = {0,
590                                        0,
591                                        0,
592                                        MatMult_Composite,
593                                        0,
594                                 /*  5*/ MatMultTranspose_Composite,
595                                        0,
596                                        0,
597                                        0,
598                                        0,
599                                 /* 10*/ 0,
600                                        0,
601                                        0,
602                                        0,
603                                        0,
604                                 /* 15*/ 0,
605                                        0,
606                                        MatGetDiagonal_Composite,
607                                        MatDiagonalScale_Composite,
608                                        0,
609                                 /* 20*/ 0,
610                                        MatAssemblyEnd_Composite,
611                                        0,
612                                        0,
613                                /* 24*/ 0,
614                                        0,
615                                        0,
616                                        0,
617                                        0,
618                                /* 29*/ 0,
619                                        0,
620                                        0,
621                                        0,
622                                        0,
623                                /* 34*/ 0,
624                                        0,
625                                        0,
626                                        0,
627                                        0,
628                                /* 39*/ 0,
629                                        0,
630                                        0,
631                                        0,
632                                        0,
633                                /* 44*/ 0,
634                                        MatScale_Composite,
635                                        MatShift_Basic,
636                                        0,
637                                        0,
638                                /* 49*/ 0,
639                                        0,
640                                        0,
641                                        0,
642                                        0,
643                                /* 54*/ 0,
644                                        0,
645                                        0,
646                                        0,
647                                        0,
648                                /* 59*/ 0,
649                                        MatDestroy_Composite,
650                                        0,
651                                        0,
652                                        0,
653                                /* 64*/ 0,
654                                        0,
655                                        0,
656                                        0,
657                                        0,
658                                /* 69*/ 0,
659                                        0,
660                                        0,
661                                        0,
662                                        0,
663                                /* 74*/ 0,
664                                        0,
665                                        0,
666                                        0,
667                                        0,
668                                /* 79*/ 0,
669                                        0,
670                                        0,
671                                        0,
672                                        0,
673                                /* 84*/ 0,
674                                        0,
675                                        0,
676                                        0,
677                                        0,
678                                /* 89*/ 0,
679                                        0,
680                                        0,
681                                        0,
682                                        0,
683                                /* 94*/ 0,
684                                        0,
685                                        0,
686                                        0,
687                                        0,
688                                 /*99*/ 0,
689                                        0,
690                                        0,
691                                        0,
692                                        0,
693                                /*104*/ 0,
694                                        0,
695                                        0,
696                                        0,
697                                        0,
698                                /*109*/ 0,
699                                        0,
700                                        0,
701                                        0,
702                                        0,
703                                /*114*/ 0,
704                                        0,
705                                        0,
706                                        0,
707                                        0,
708                                /*119*/ 0,
709                                        0,
710                                        0,
711                                        0,
712                                        0,
713                                /*124*/ 0,
714                                        0,
715                                        0,
716                                        0,
717                                        0,
718                                /*129*/ 0,
719                                        0,
720                                        0,
721                                        0,
722                                        0,
723                                /*134*/ 0,
724                                        0,
725                                        0,
726                                        0,
727                                        0,
728                                /*139*/ 0,
729                                        0,
730                                        0
731 };
732 
733 /*MC
734    MATCOMPOSITE - A matrix defined by the sum (or product) of one or more matrices.
735     The matrices need to have a correct size and parallel layout for the sum or product to be valid.
736 
737    Notes:
738     to use the product of the matrices call MatCompositeSetType(mat,MAT_COMPOSITE_MULTIPLICATIVE);
739 
740   Level: advanced
741 
742 .seealso: MatCreateComposite(), MatCompositeAddMat(), MatSetType(), MatCompositeSetType(), MatCompositeMerge(), MatCompositeSetMergeFromRight(), MatCompositeGetNmat(), MatCompositeGetMat()
743 M*/
744 
745 PETSC_EXTERN PetscErrorCode MatCreate_Composite(Mat A)
746 {
747   Mat_Composite  *b;
748   PetscErrorCode ierr;
749 
750   PetscFunctionBegin;
751   ierr    = PetscNewLog(A,&b);CHKERRQ(ierr);
752   A->data = (void*)b;
753   ierr    = PetscMemcpy(A->ops,&MatOps_Values,sizeof(struct _MatOps));CHKERRQ(ierr);
754 
755   ierr = PetscLayoutSetUp(A->rmap);CHKERRQ(ierr);
756   ierr = PetscLayoutSetUp(A->cmap);CHKERRQ(ierr);
757 
758   A->assembled      = PETSC_TRUE;
759   A->preallocated   = PETSC_TRUE;
760   b->type           = MAT_COMPOSITE_ADDITIVE;
761   b->scale          = 1.0;
762   b->nmat           = 0;
763   b->mergefromright = PETSC_TRUE;
764 
765   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeAddMat_C",MatCompositeAddMat_Composite);CHKERRQ(ierr);
766   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeSetType_C",MatCompositeSetType_Composite);CHKERRQ(ierr);
767   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeMerge_C",MatCompositeMerge_Composite);CHKERRQ(ierr);
768   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeGetNMat_C",MatCompositeGetNMat_Composite);CHKERRQ(ierr);
769   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeGetMat_C",MatCompositeGetMat_Composite);CHKERRQ(ierr);
770   ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeSetMergeFromRight_C",MatCompositeSetMergeFromRight_Composite);CHKERRQ(ierr);
771 
772   ierr = PetscObjectChangeTypeName((PetscObject)A,MATCOMPOSITE);CHKERRQ(ierr);
773   PetscFunctionReturn(0);
774 }
775 
776