Actual source code: pack.c

petsc-3.6.1 2015-08-06
Report Typos and Errors
  2: #include <../src/dm/impls/composite/packimpl.h>       /*I  "petscdmcomposite.h"  I*/
  3: #include <petsc/private/isimpl.h>
  4: #include <petscds.h>

  8: /*@C
  9:     DMCompositeSetCoupling - Sets user provided routines that compute the coupling between the
 10:       separate components (DMs) in a DMto build the correct matrix nonzero structure.


 13:     Logically Collective on MPI_Comm

 15:     Input Parameter:
 16: +   dm - the composite object
 17: -   formcouplelocations - routine to set the nonzero locations in the matrix

 19:     Level: advanced

 21:     Notes: See DMSetApplicationContext() and DMGetApplicationContext() for how to get user information into
 22:         this routine

 24: @*/
 25: PetscErrorCode  DMCompositeSetCoupling(DM dm,PetscErrorCode (*FormCoupleLocations)(DM,Mat,PetscInt*,PetscInt*,PetscInt,PetscInt,PetscInt,PetscInt))
 26: {
 27:   DM_Composite *com = (DM_Composite*)dm->data;

 30:   com->FormCoupleLocations = FormCoupleLocations;
 31:   return(0);
 32: }

 36: PetscErrorCode  DMDestroy_Composite(DM dm)
 37: {
 38:   PetscErrorCode         ierr;
 39:   struct DMCompositeLink *next, *prev;
 40:   DM_Composite           *com = (DM_Composite*)dm->data;

 43:   next = com->next;
 44:   while (next) {
 45:     prev = next;
 46:     next = next->next;
 47:     DMDestroy(&prev->dm);
 48:     PetscFree(prev->grstarts);
 49:     PetscFree(prev);
 50:   }
 51:   /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
 52:   PetscFree(com);
 53:   return(0);
 54: }

 58: PetscErrorCode  DMView_Composite(DM dm,PetscViewer v)
 59: {
 61:   PetscBool      iascii;
 62:   DM_Composite   *com = (DM_Composite*)dm->data;

 65:   PetscObjectTypeCompare((PetscObject)v,PETSCVIEWERASCII,&iascii);
 66:   if (iascii) {
 67:     struct DMCompositeLink *lnk = com->next;
 68:     PetscInt               i;

 70:     PetscViewerASCIIPrintf(v,"DM (%s)\n",((PetscObject)dm)->prefix ? ((PetscObject)dm)->prefix : "no prefix");
 71:     PetscViewerASCIIPrintf(v,"  contains %D DMs\n",com->nDM);
 72:     PetscViewerASCIIPushTab(v);
 73:     for (i=0; lnk; lnk=lnk->next,i++) {
 74:       PetscViewerASCIIPrintf(v,"Link %D: DM of type %s\n",i,((PetscObject)lnk->dm)->type_name);
 75:       PetscViewerASCIIPushTab(v);
 76:       DMView(lnk->dm,v);
 77:       PetscViewerASCIIPopTab(v);
 78:     }
 79:     PetscViewerASCIIPopTab(v);
 80:   }
 81:   return(0);
 82: }

 84: /* --------------------------------------------------------------------------------------*/
 87: PetscErrorCode  DMSetUp_Composite(DM dm)
 88: {
 89:   PetscErrorCode         ierr;
 90:   PetscInt               nprev = 0;
 91:   PetscMPIInt            rank,size;
 92:   DM_Composite           *com  = (DM_Composite*)dm->data;
 93:   struct DMCompositeLink *next = com->next;
 94:   PetscLayout            map;

 97:   if (com->setup) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Packer has already been setup");
 98:   PetscLayoutCreate(PetscObjectComm((PetscObject)dm),&map);
 99:   PetscLayoutSetLocalSize(map,com->n);
100:   PetscLayoutSetSize(map,PETSC_DETERMINE);
101:   PetscLayoutSetBlockSize(map,1);
102:   PetscLayoutSetUp(map);
103:   PetscLayoutGetSize(map,&com->N);
104:   PetscLayoutGetRange(map,&com->rstart,NULL);
105:   PetscLayoutDestroy(&map);

107:   /* now set the rstart for each linked vector */
108:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);
109:   MPI_Comm_size(PetscObjectComm((PetscObject)dm),&size);
110:   while (next) {
111:     next->rstart  = nprev;
112:     nprev        += next->n;
113:     next->grstart = com->rstart + next->rstart;
114:     PetscMalloc1(size,&next->grstarts);
115:     MPI_Allgather(&next->grstart,1,MPIU_INT,next->grstarts,1,MPIU_INT,PetscObjectComm((PetscObject)dm));
116:     next          = next->next;
117:   }
118:   com->setup = PETSC_TRUE;
119:   return(0);
120: }

122: /* ----------------------------------------------------------------------------------*/

126: /*@
127:     DMCompositeGetNumberDM - Get's the number of DM objects in the DMComposite
128:        representation.

130:     Not Collective

132:     Input Parameter:
133: .    dm - the packer object

135:     Output Parameter:
136: .     nDM - the number of DMs

138:     Level: beginner

140: @*/
141: PetscErrorCode  DMCompositeGetNumberDM(DM dm,PetscInt *nDM)
142: {
143:   DM_Composite *com = (DM_Composite*)dm->data;

147:   *nDM = com->nDM;
148:   return(0);
149: }


154: /*@C
155:     DMCompositeGetAccess - Allows one to access the individual packed vectors in their global
156:        representation.

158:     Collective on DMComposite

160:     Input Parameters:
161: +    dm - the packer object
162: -    gvec - the global vector

164:     Output Parameters:
165: .    Vec* ... - the packed parallel vectors, NULL for those that are not needed

167:     Notes: Use DMCompositeRestoreAccess() to return the vectors when you no longer need them

169:     Fortran Notes:

171:     Fortran callers must use numbered versions of this routine, e.g., DMCompositeGetAccess4(dm,gvec,vec1,vec2,vec3,vec4)
172:     or use the alternative interface DMCompositeGetAccessArray().

174:     Level: advanced

176: .seealso: DMCompositeGetEntries(), DMCompositeScatter()
177: @*/
178: PetscErrorCode  DMCompositeGetAccess(DM dm,Vec gvec,...)
179: {
180:   va_list                Argp;
181:   PetscErrorCode         ierr;
182:   struct DMCompositeLink *next;
183:   DM_Composite           *com = (DM_Composite*)dm->data;
184:   PetscInt               readonly;

189:   next = com->next;
190:   if (!com->setup) {
191:     DMSetUp(dm);
192:   }

194:   VecLockGet(gvec,&readonly);
195:   /* loop over packed objects, handling one at at time */
196:   va_start(Argp,gvec);
197:   while (next) {
198:     Vec *vec;
199:     vec = va_arg(Argp, Vec*);
200:     if (vec) {
201:       DMGetGlobalVector(next->dm,vec);
202:       if (readonly) {
203:         const PetscScalar *array;
204:         VecGetArrayRead(gvec,&array);
205:         VecPlaceArray(*vec,array+next->rstart);
206:         VecLockPush(*vec);
207:         VecRestoreArrayRead(gvec,&array);
208:       } else {
209:         PetscScalar *array;
210:         VecGetArray(gvec,&array);
211:         VecPlaceArray(*vec,array+next->rstart);
212:         VecRestoreArray(gvec,&array);
213:       }
214:     }
215:     next = next->next;
216:   }
217:   va_end(Argp);
218:   return(0);
219: }

223: /*@C
224:     DMCompositeGetAccessArray - Allows one to access the individual packed vectors in their global
225:        representation.

227:     Collective on DMComposite

229:     Input Parameters:
230: +    dm - the packer object
231: .    pvec - packed vector
232: .    nwanted - number of vectors wanted
233: -    wanted - sorted array of vectors wanted, or NULL to get all vectors

235:     Output Parameters:
236: .    vecs - array of requested global vectors (must be allocated)

238:     Notes: Use DMCompositeRestoreAccessArray() to return the vectors when you no longer need them

240:     Level: advanced

242: .seealso: DMCompositeGetAccess(), DMCompositeGetEntries(), DMCompositeScatter(), DMCompositeGather()
243: @*/
244: PetscErrorCode  DMCompositeGetAccessArray(DM dm,Vec pvec,PetscInt nwanted,const PetscInt *wanted,Vec *vecs)
245: {
246:   PetscErrorCode         ierr;
247:   struct DMCompositeLink *link;
248:   PetscInt               i,wnum;
249:   DM_Composite           *com = (DM_Composite*)dm->data;
250:   PetscInt               readonly;
251: 
255:   if (!com->setup) {
256:     DMSetUp(dm);
257:   }

259:   VecLockGet(pvec,&readonly);
260:   for (i=0,wnum=0,link=com->next; link && wnum<nwanted; i++,link=link->next) {
261:     if (!wanted || i == wanted[wnum]) {
262:       Vec v;
263:       DMGetGlobalVector(link->dm,&v);
264:       if (readonly) {
265:         const PetscScalar *array;
266:         VecGetArrayRead(pvec,&array);
267:         VecPlaceArray(v,array+link->rstart);
268:         VecLockPush(v);
269:         VecRestoreArrayRead(pvec,&array);
270:       } else {
271:         PetscScalar *array;
272:         VecGetArray(pvec,&array);
273:         VecPlaceArray(v,array+link->rstart);
274:         VecRestoreArray(pvec,&array);
275:       }
276:       vecs[wnum++] = v;
277:     }
278:   }
279:   return(0);
280: }

284: /*@C
285:     DMCompositeRestoreAccess - Returns the vectors obtained with DMCompositeGetAccess()
286:        representation.

288:     Collective on DMComposite

290:     Input Parameters:
291: +    dm - the packer object
292: .    gvec - the global vector
293: -    Vec* ... - the individual parallel vectors, NULL for those that are not needed

295:     Level: advanced

297: .seealso  DMCompositeAddDM(), DMCreateGlobalVector(),
298:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeScatter(),
299:          DMCompositeRestoreAccess(), DMCompositeGetAccess()

301: @*/
302: PetscErrorCode  DMCompositeRestoreAccess(DM dm,Vec gvec,...)
303: {
304:   va_list                Argp;
305:   PetscErrorCode         ierr;
306:   struct DMCompositeLink *next;
307:   DM_Composite           *com = (DM_Composite*)dm->data;
308:   PetscInt               readonly;

313:   next = com->next;
314:   if (!com->setup) {
315:     DMSetUp(dm);
316:   }

318:   VecLockGet(gvec,&readonly);
319:   /* loop over packed objects, handling one at at time */
320:   va_start(Argp,gvec);
321:   while (next) {
322:     Vec *vec;
323:     vec = va_arg(Argp, Vec*);
324:     if (vec) {
325:       VecResetArray(*vec);
326:       if (readonly) {
327:         VecLockPop(*vec);
328:       }
329:       DMRestoreGlobalVector(next->dm,vec);
330:     }
331:     next = next->next;
332:   }
333:   va_end(Argp);
334:   return(0);
335: }

339: /*@C
340:     DMCompositeRestoreAccessArray - Returns the vectors obtained with DMCompositeGetAccessArray()

342:     Collective on DMComposite

344:     Input Parameters:
345: +    dm - the packer object
346: .    pvec - packed vector
347: .    nwanted - number of vectors wanted
348: .    wanted - sorted array of vectors wanted, or NULL to get all vectors
349: -    vecs - array of global vectors to return

351:     Level: advanced

353: .seealso: DMCompositeRestoreAccess(), DMCompositeRestoreEntries(), DMCompositeScatter(), DMCompositeGather()
354: @*/
355: PetscErrorCode  DMCompositeRestoreAccessArray(DM dm,Vec pvec,PetscInt nwanted,const PetscInt *wanted,Vec *vecs)
356: {
357:   PetscErrorCode         ierr;
358:   struct DMCompositeLink *link;
359:   PetscInt               i,wnum;
360:   DM_Composite           *com = (DM_Composite*)dm->data;
361:   PetscInt               readonly;

366:   if (!com->setup) {
367:     DMSetUp(dm);
368:   }

370:   VecLockGet(pvec,&readonly);
371:   for (i=0,wnum=0,link=com->next; link && wnum<nwanted; i++,link=link->next) {
372:     if (!wanted || i == wanted[wnum]) {
373:       VecResetArray(vecs[wnum]);
374:       if (readonly) {
375:         VecLockPop(vecs[wnum]);
376:       }
377:       DMRestoreGlobalVector(link->dm,&vecs[wnum]);
378:       wnum++;
379:     }
380:   }
381:   return(0);
382: }

386: /*@C
387:     DMCompositeScatter - Scatters from a global packed vector into its individual local vectors

389:     Collective on DMComposite

391:     Input Parameters:
392: +    dm - the packer object
393: .    gvec - the global vector
394: -    Vec ... - the individual sequential vectors, NULL for those that are not needed

396:     Level: advanced

398:     Notes:
399:     DMCompositeScatterArray() is a non-variadic alternative that is often more convenient for library callers and is
400:     accessible from Fortran.

402: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
403:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
404:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries()
405:          DMCompositeScatterArray()

407: @*/
408: PetscErrorCode  DMCompositeScatter(DM dm,Vec gvec,...)
409: {
410:   va_list                Argp;
411:   PetscErrorCode         ierr;
412:   struct DMCompositeLink *next;
413:   PetscInt               cnt;
414:   DM_Composite           *com = (DM_Composite*)dm->data;

419:   if (!com->setup) {
420:     DMSetUp(dm);
421:   }

423:   /* loop over packed objects, handling one at at time */
424:   va_start(Argp,gvec);
425:   for (cnt=3,next=com->next; next; cnt++,next=next->next) {
426:     Vec local;
427:     local = va_arg(Argp, Vec);
428:     if (local) {
429:       Vec               global;
430:       const PetscScalar *array;
432:       DMGetGlobalVector(next->dm,&global);
433:       VecGetArrayRead(gvec,&array);
434:       VecPlaceArray(global,array+next->rstart);
435:       DMGlobalToLocalBegin(next->dm,global,INSERT_VALUES,local);
436:       DMGlobalToLocalEnd(next->dm,global,INSERT_VALUES,local);
437:       VecRestoreArrayRead(gvec,&array);
438:       VecResetArray(global);
439:       DMRestoreGlobalVector(next->dm,&global);
440:     }
441:   }
442:   va_end(Argp);
443:   return(0);
444: }

448: /*@
449:     DMCompositeScatterArray - Scatters from a global packed vector into its individual local vectors

451:     Collective on DMComposite

453:     Input Parameters:
454: +    dm - the packer object
455: .    gvec - the global vector
456: .    lvecs - array of local vectors, NULL for any that are not needed

458:     Level: advanced

460:     Note:
461:     This is a non-variadic alternative to DMCompositeScatter()

463: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector()
464:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
465:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries()

467: @*/
468: PetscErrorCode  DMCompositeScatterArray(DM dm,Vec gvec,Vec *lvecs)
469: {
470:   PetscErrorCode         ierr;
471:   struct DMCompositeLink *next;
472:   PetscInt               i;
473:   DM_Composite           *com = (DM_Composite*)dm->data;

478:   if (!com->setup) {
479:     DMSetUp(dm);
480:   }

482:   /* loop over packed objects, handling one at at time */
483:   for (i=0,next=com->next; next; next=next->next,i++) {
484:     if (lvecs[i]) {
485:       Vec         global;
486:       const PetscScalar *array;
488:       DMGetGlobalVector(next->dm,&global);
489:       VecGetArrayRead(gvec,&array);
490:       VecPlaceArray(global,(PetscScalar*)array+next->rstart);
491:       DMGlobalToLocalBegin(next->dm,global,INSERT_VALUES,lvecs[i]);
492:       DMGlobalToLocalEnd(next->dm,global,INSERT_VALUES,lvecs[i]);
493:       VecRestoreArrayRead(gvec,&array);
494:       VecResetArray(global);
495:       DMRestoreGlobalVector(next->dm,&global);
496:     }
497:   }
498:   return(0);
499: }

503: /*@C
504:     DMCompositeGather - Gathers into a global packed vector from its individual local vectors

506:     Collective on DMComposite

508:     Input Parameter:
509: +    dm - the packer object
510: .    gvec - the global vector
511: .    imode - INSERT_VALUES or ADD_VALUES
512: -    Vec ... - the individual sequential vectors, NULL for any that are not needed

514:     Level: advanced

516: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
517:          DMCompositeScatter(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
518:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries()

520: @*/
521: PetscErrorCode  DMCompositeGather(DM dm,Vec gvec,InsertMode imode,...)
522: {
523:   va_list                Argp;
524:   PetscErrorCode         ierr;
525:   struct DMCompositeLink *next;
526:   DM_Composite           *com = (DM_Composite*)dm->data;
527:   PetscInt               cnt;

532:   if (!com->setup) {
533:     DMSetUp(dm);
534:   }

536:   /* loop over packed objects, handling one at at time */
537:   va_start(Argp,imode);
538:   for (cnt=3,next=com->next; next; cnt++,next=next->next) {
539:     Vec local;
540:     local = va_arg(Argp, Vec);
541:     if (local) {
542:       PetscScalar *array;
543:       Vec         global;
545:       DMGetGlobalVector(next->dm,&global);
546:       VecGetArray(gvec,&array);
547:       VecPlaceArray(global,array+next->rstart);
548:       DMLocalToGlobalBegin(next->dm,local,imode,global);
549:       DMLocalToGlobalEnd(next->dm,local,imode,global);
550:       VecRestoreArray(gvec,&array);
551:       VecResetArray(global);
552:       DMRestoreGlobalVector(next->dm,&global);
553:     }
554:   }
555:   va_end(Argp);
556:   return(0);
557: }

561: /*@
562:     DMCompositeGatherArray - Gathers into a global packed vector from its individual local vectors

564:     Collective on DMComposite

566:     Input Parameter:
567: +    dm - the packer object
568: .    gvec - the global vector
569: .    imode - INSERT_VALUES or ADD_VALUES
570: -    lvecs - the individual sequential vectors, NULL for any that are not needed

572:     Level: advanced

574:     Notes:
575:     This is a non-variadic alternative to DMCompositeGather().

577: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
578:          DMCompositeScatter(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
579:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries(),
580: @*/
581: PetscErrorCode  DMCompositeGatherArray(DM dm,Vec gvec,InsertMode imode,Vec *lvecs)
582: {
583:   PetscErrorCode         ierr;
584:   struct DMCompositeLink *next;
585:   DM_Composite           *com = (DM_Composite*)dm->data;
586:   PetscInt               i;

591:   if (!com->setup) {
592:     DMSetUp(dm);
593:   }

595:   /* loop over packed objects, handling one at at time */
596:   for (next=com->next,i=0; next; next=next->next,i++) {
597:     if (lvecs[i]) {
598:       PetscScalar *array;
599:       Vec         global;
601:       DMGetGlobalVector(next->dm,&global);
602:       VecGetArray(gvec,&array);
603:       VecPlaceArray(global,array+next->rstart);
604:       DMLocalToGlobalBegin(next->dm,lvecs[i],imode,global);
605:       DMLocalToGlobalEnd(next->dm,lvecs[i],imode,global);
606:       VecRestoreArray(gvec,&array);
607:       VecResetArray(global);
608:       DMRestoreGlobalVector(next->dm,&global);
609:     }
610:   }
611:   return(0);
612: }

616: /*@C
617:     DMCompositeAddDM - adds a DM  vector to a DMComposite

619:     Collective on DMComposite

621:     Input Parameter:
622: +    dm - the packer object
623: -    dm - the DM object, if the DM is a da you will need to caste it with a (DM)

625:     Level: advanced

627: .seealso DMDestroy(), DMCompositeGather(), DMCompositeAddDM(), DMCreateGlobalVector(),
628:          DMCompositeScatter(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
629:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries()

631: @*/
632: PetscErrorCode  DMCompositeAddDM(DM dmc,DM dm)
633: {
634:   PetscErrorCode         ierr;
635:   PetscInt               n,nlocal;
636:   struct DMCompositeLink *mine,*next;
637:   Vec                    global,local;
638:   DM_Composite           *com = (DM_Composite*)dmc->data;

643:   next = com->next;
644:   if (com->setup) SETERRQ(PetscObjectComm((PetscObject)dmc),PETSC_ERR_ARG_WRONGSTATE,"Cannot add a DM once you have used the DMComposite");

646:   /* create new link */
647:   PetscNew(&mine);
648:   PetscObjectReference((PetscObject)dm);
649:   DMGetGlobalVector(dm,&global);
650:   VecGetLocalSize(global,&n);
651:   DMRestoreGlobalVector(dm,&global);
652:   DMGetLocalVector(dm,&local);
653:   VecGetSize(local,&nlocal);
654:   DMRestoreLocalVector(dm,&local);

656:   mine->n      = n;
657:   mine->nlocal = nlocal;
658:   mine->dm     = dm;
659:   mine->next   = NULL;
660:   com->n      += n;

662:   /* add to end of list */
663:   if (!next) com->next = mine;
664:   else {
665:     while (next->next) next = next->next;
666:     next->next = mine;
667:   }
668:   com->nDM++;
669:   com->nmine++;
670:   return(0);
671: }

673: #include <petscdraw.h>
674: PETSC_EXTERN PetscErrorCode  VecView_MPI(Vec,PetscViewer);
677: PetscErrorCode  VecView_DMComposite(Vec gvec,PetscViewer viewer)
678: {
679:   DM                     dm;
680:   PetscErrorCode         ierr;
681:   struct DMCompositeLink *next;
682:   PetscBool              isdraw;
683:   DM_Composite           *com;

686:   VecGetDM(gvec, &dm);
687:   if (!dm) SETERRQ(PetscObjectComm((PetscObject)gvec),PETSC_ERR_ARG_WRONG,"Vector not generated from a DMComposite");
688:   com  = (DM_Composite*)dm->data;
689:   next = com->next;

691:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
692:   if (!isdraw) {
693:     /* do I really want to call this? */
694:     VecView_MPI(gvec,viewer);
695:   } else {
696:     PetscInt cnt = 0;

698:     /* loop over packed objects, handling one at at time */
699:     while (next) {
700:       Vec         vec;
701:       PetscScalar *array;
702:       PetscInt    bs;

704:       /* Should use VecGetSubVector() eventually, but would need to forward the DM for that to work */
705:       DMGetGlobalVector(next->dm,&vec);
706:       VecGetArray(gvec,&array);
707:       VecPlaceArray(vec,array+next->rstart);
708:       VecRestoreArray(gvec,&array);
709:       VecView(vec,viewer);
710:       VecGetBlockSize(vec,&bs);
711:       VecResetArray(vec);
712:       DMRestoreGlobalVector(next->dm,&vec);
713:       PetscViewerDrawBaseAdd(viewer,bs);
714:       cnt += bs;
715:       next = next->next;
716:     }
717:     PetscViewerDrawBaseAdd(viewer,-cnt);
718:   }
719:   return(0);
720: }

724: PetscErrorCode  DMCreateGlobalVector_Composite(DM dm,Vec *gvec)
725: {
727:   DM_Composite   *com = (DM_Composite*)dm->data;

731:   DMSetUp(dm);
732:   VecCreateMPI(PetscObjectComm((PetscObject)dm),com->n,com->N,gvec);
733:   VecSetDM(*gvec, dm);
734:   VecSetOperation(*gvec,VECOP_VIEW,(void (*)(void))VecView_DMComposite);
735:   return(0);
736: }

740: PetscErrorCode  DMCreateLocalVector_Composite(DM dm,Vec *lvec)
741: {
743:   DM_Composite   *com = (DM_Composite*)dm->data;

747:   if (!com->setup) {
748:     DMSetUp(dm);
749:   }
750:   VecCreateSeq(PetscObjectComm((PetscObject)dm),com->nghost,lvec);
751:   VecSetDM(*lvec, dm);
752:   return(0);
753: }

757: /*@C
758:     DMCompositeGetISLocalToGlobalMappings - gets an ISLocalToGlobalMapping for each DM in the DMComposite, maps to the composite global space

760:     Collective on DM

762:     Input Parameter:
763: .    dm - the packer object

765:     Output Parameters:
766: .    ltogs - the individual mappings for each packed vector. Note that this includes
767:            all the ghost points that individual ghosted DMDA's may have.

769:     Level: advanced

771:     Notes:
772:        Each entry of ltogs should be destroyed with ISLocalToGlobalMappingDestroy(), the ltogs array should be freed with PetscFree().

774: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
775:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetAccess(), DMCompositeScatter(),
776:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(),DMCompositeGetEntries()

778: @*/
779: PetscErrorCode  DMCompositeGetISLocalToGlobalMappings(DM dm,ISLocalToGlobalMapping **ltogs)
780: {
781:   PetscErrorCode         ierr;
782:   PetscInt               i,*idx,n,cnt;
783:   struct DMCompositeLink *next;
784:   PetscMPIInt            rank;
785:   DM_Composite           *com = (DM_Composite*)dm->data;

789:   DMSetUp(dm);
790:   PetscMalloc1(com->nDM,ltogs);
791:   next = com->next;
792:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);

794:   /* loop over packed objects, handling one at at time */
795:   cnt = 0;
796:   while (next) {
797:     ISLocalToGlobalMapping ltog;
798:     PetscMPIInt            size;
799:     const PetscInt         *suboff,*indices;
800:     Vec                    global;

802:     /* Get sub-DM global indices for each local dof */
803:     DMGetLocalToGlobalMapping(next->dm,&ltog);
804:     ISLocalToGlobalMappingGetSize(ltog,&n);
805:     ISLocalToGlobalMappingGetIndices(ltog,&indices);
806:     PetscMalloc1(n,&idx);

808:     /* Get the offsets for the sub-DM global vector */
809:     DMGetGlobalVector(next->dm,&global);
810:     VecGetOwnershipRanges(global,&suboff);
811:     MPI_Comm_size(PetscObjectComm((PetscObject)global),&size);

813:     /* Shift the sub-DM definition of the global space to the composite global space */
814:     for (i=0; i<n; i++) {
815:       PetscInt subi = indices[i],lo = 0,hi = size,t;
816:       /* Binary search to find which rank owns subi */
817:       while (hi-lo > 1) {
818:         t = lo + (hi-lo)/2;
819:         if (suboff[t] > subi) hi = t;
820:         else                  lo = t;
821:       }
822:       idx[i] = subi - suboff[lo] + next->grstarts[lo];
823:     }
824:     ISLocalToGlobalMappingRestoreIndices(ltog,&indices);
825:     ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)dm),1,n,idx,PETSC_OWN_POINTER,&(*ltogs)[cnt]);
826:     DMRestoreGlobalVector(next->dm,&global);
827:     next = next->next;
828:     cnt++;
829:   }
830:   return(0);
831: }

835: /*@C
836:    DMCompositeGetLocalISs - Gets index sets for each component of a composite local vector

838:    Not Collective

840:    Input Arguments:
841: . dm - composite DM

843:    Output Arguments:
844: . is - array of serial index sets for each each component of the DMComposite

846:    Level: intermediate

848:    Notes:
849:    At present, a composite local vector does not normally exist.  This function is used to provide index sets for
850:    MatGetLocalSubMatrix().  In the future, the scatters for each entry in the DMComposite may be be merged into a single
851:    scatter to a composite local vector.  The user should not typically need to know which is being done.

853:    To get the composite global indices at all local points (including ghosts), use DMCompositeGetISLocalToGlobalMappings().

855:    To get index sets for pieces of the composite global vector, use DMCompositeGetGlobalISs().

857:    Each returned IS should be destroyed with ISDestroy(), the array should be freed with PetscFree().

859: .seealso: DMCompositeGetGlobalISs(), DMCompositeGetISLocalToGlobalMappings(), MatGetLocalSubMatrix(), MatCreateLocalRef()
860: @*/
861: PetscErrorCode  DMCompositeGetLocalISs(DM dm,IS **is)
862: {
863:   PetscErrorCode         ierr;
864:   DM_Composite           *com = (DM_Composite*)dm->data;
865:   struct DMCompositeLink *link;
866:   PetscInt               cnt,start;

871:   PetscMalloc1(com->nmine,is);
872:   for (cnt=0,start=0,link=com->next; link; start+=link->nlocal,cnt++,link=link->next) {
873:     PetscInt bs;
874:     ISCreateStride(PETSC_COMM_SELF,link->nlocal,start,1,&(*is)[cnt]);
875:     DMGetBlockSize(link->dm,&bs);
876:     ISSetBlockSize((*is)[cnt],bs);
877:   }
878:   return(0);
879: }

883: /*@C
884:     DMCompositeGetGlobalISs - Gets the index sets for each composed object

886:     Collective on DMComposite

888:     Input Parameter:
889: .    dm - the packer object

891:     Output Parameters:
892: .    is - the array of index sets

894:     Level: advanced

896:     Notes:
897:        The is entries should be destroyed with ISDestroy(), the is array should be freed with PetscFree()

899:        These could be used to extract a subset of vector entries for a "multi-physics" preconditioner

901:        Use DMCompositeGetLocalISs() for index sets in the packed local numbering, and
902:        DMCompositeGetISLocalToGlobalMappings() for to map local sub-DM (including ghost) indices to packed global
903:        indices.

905: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
906:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetAccess(), DMCompositeScatter(),
907:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(),DMCompositeGetEntries()

909: @*/

911: PetscErrorCode  DMCompositeGetGlobalISs(DM dm,IS *is[])
912: {
913:   PetscErrorCode         ierr;
914:   PetscInt               cnt = 0;
915:   struct DMCompositeLink *next;
916:   PetscMPIInt            rank;
917:   DM_Composite           *com = (DM_Composite*)dm->data;

921:   PetscMalloc1(com->nDM,is);
922:   next = com->next;
923:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);

925:   /* loop over packed objects, handling one at at time */
926:   while (next) {
927:     ISCreateStride(PetscObjectComm((PetscObject)dm),next->n,next->grstart,1,&(*is)[cnt]);
928:     if (dm->prob) {
929:       MatNullSpace space;
930:       Mat          pmat;
931:       PetscObject  disc;
932:       PetscInt     Nf;

934:       PetscDSGetNumFields(dm->prob, &Nf);
935:       if (cnt < Nf) {
936:         PetscDSGetDiscretization(dm->prob, cnt, &disc);
937:         PetscObjectQuery(disc, "nullspace", (PetscObject*) &space);
938:         if (space) {PetscObjectCompose((PetscObject) (*is)[cnt], "nullspace", (PetscObject) space);}
939:         PetscObjectQuery(disc, "nearnullspace", (PetscObject*) &space);
940:         if (space) {PetscObjectCompose((PetscObject) (*is)[cnt], "nearnullspace", (PetscObject) space);}
941:         PetscObjectQuery(disc, "pmat", (PetscObject*) &pmat);
942:         if (pmat) {PetscObjectCompose((PetscObject) (*is)[cnt], "pmat", (PetscObject) pmat);}
943:       }
944:     }
945:     cnt++;
946:     next = next->next;
947:   }
948:   return(0);
949: }

953: PetscErrorCode DMCreateFieldIS_Composite(DM dm, PetscInt *numFields,char ***fieldNames, IS **fields)
954: {
955:   PetscInt       nDM;
956:   DM             *dms;
957:   PetscInt       i;

961:   DMCompositeGetNumberDM(dm, &nDM);
962:   if (numFields) *numFields = nDM;
963:   DMCompositeGetGlobalISs(dm, fields);
964:   if (fieldNames) {
965:     PetscMalloc1(nDM, &dms);
966:     PetscMalloc1(nDM, fieldNames);
967:     DMCompositeGetEntriesArray(dm, dms);
968:     for (i=0; i<nDM; i++) {
969:       char       buf[256];
970:       const char *splitname;

972:       /* Split naming precedence: object name, prefix, number */
973:       splitname = ((PetscObject) dm)->name;
974:       if (!splitname) {
975:         PetscObjectGetOptionsPrefix((PetscObject)dms[i],&splitname);
976:         if (splitname) {
977:           size_t len;
978:           PetscStrncpy(buf,splitname,sizeof(buf));
979:           buf[sizeof(buf) - 1] = 0;
980:           PetscStrlen(buf,&len);
981:           if (buf[len-1] == '_') buf[len-1] = 0; /* Remove trailing underscore if it was used */
982:           splitname = buf;
983:         }
984:       }
985:       if (!splitname) {
986:         PetscSNPrintf(buf,sizeof(buf),"%D",i);
987:         splitname = buf;
988:       }
989:       PetscStrallocpy(splitname,&(*fieldNames)[i]);
990:     }
991:     PetscFree(dms);
992:   }
993:   return(0);
994: }

996: /*
997:  This could take over from DMCreateFieldIS(), as it is more general,
998:  making DMCreateFieldIS() a special case -- calling with dmlist == NULL;
999:  At this point it's probably best to be less intrusive, however.
1000:  */
1003: PetscErrorCode DMCreateFieldDecomposition_Composite(DM dm, PetscInt *len,char ***namelist, IS **islist, DM **dmlist)
1004: {
1005:   PetscInt       nDM;
1006:   PetscInt       i;

1010:   DMCreateFieldIS_Composite(dm, len, namelist, islist);
1011:   if (dmlist) {
1012:     DMCompositeGetNumberDM(dm, &nDM);
1013:     PetscMalloc1(nDM, dmlist);
1014:     DMCompositeGetEntriesArray(dm, *dmlist);
1015:     for (i=0; i<nDM; i++) {
1016:       PetscObjectReference((PetscObject)((*dmlist)[i]));
1017:     }
1018:   }
1019:   return(0);
1020: }



1024: /* -------------------------------------------------------------------------------------*/
1027: /*@C
1028:     DMCompositeGetLocalVectors - Gets local vectors for each part of a DMComposite.
1029:        Use DMCompositeRestoreLocalVectors() to return them.

1031:     Not Collective

1033:     Input Parameter:
1034: .    dm - the packer object

1036:     Output Parameter:
1037: .   Vec ... - the individual sequential Vecs

1039:     Level: advanced

1041: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
1042:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
1043:          DMCompositeRestoreLocalVectors(), DMCompositeScatter(), DMCompositeGetEntries()

1045: @*/
1046: PetscErrorCode  DMCompositeGetLocalVectors(DM dm,...)
1047: {
1048:   va_list                Argp;
1049:   PetscErrorCode         ierr;
1050:   struct DMCompositeLink *next;
1051:   DM_Composite           *com = (DM_Composite*)dm->data;

1055:   next = com->next;
1056:   /* loop over packed objects, handling one at at time */
1057:   va_start(Argp,dm);
1058:   while (next) {
1059:     Vec *vec;
1060:     vec = va_arg(Argp, Vec*);
1061:     if (vec) {DMGetLocalVector(next->dm,vec);}
1062:     next = next->next;
1063:   }
1064:   va_end(Argp);
1065:   return(0);
1066: }

1070: /*@C
1071:     DMCompositeRestoreLocalVectors - Restores local vectors for each part of a DMComposite.

1073:     Not Collective

1075:     Input Parameter:
1076: .    dm - the packer object

1078:     Output Parameter:
1079: .   Vec ... - the individual sequential Vecs

1081:     Level: advanced

1083: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(),
1084:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
1085:          DMCompositeGetLocalVectors(), DMCompositeScatter(), DMCompositeGetEntries()

1087: @*/
1088: PetscErrorCode  DMCompositeRestoreLocalVectors(DM dm,...)
1089: {
1090:   va_list                Argp;
1091:   PetscErrorCode         ierr;
1092:   struct DMCompositeLink *next;
1093:   DM_Composite           *com = (DM_Composite*)dm->data;

1097:   next = com->next;
1098:   /* loop over packed objects, handling one at at time */
1099:   va_start(Argp,dm);
1100:   while (next) {
1101:     Vec *vec;
1102:     vec = va_arg(Argp, Vec*);
1103:     if (vec) {DMRestoreLocalVector(next->dm,vec);}
1104:     next = next->next;
1105:   }
1106:   va_end(Argp);
1107:   return(0);
1108: }

1110: /* -------------------------------------------------------------------------------------*/
1113: /*@C
1114:     DMCompositeGetEntries - Gets the DM for each entry in a DMComposite.

1116:     Not Collective

1118:     Input Parameter:
1119: .    dm - the packer object

1121:     Output Parameter:
1122: .   DM ... - the individual entries (DMs)

1124:     Level: advanced

1126: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(), DMCompositeGetEntriesArray()
1127:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
1128:          DMCompositeRestoreLocalVectors(), DMCompositeGetLocalVectors(),  DMCompositeScatter(),
1129:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors()

1131: @*/
1132: PetscErrorCode  DMCompositeGetEntries(DM dm,...)
1133: {
1134:   va_list                Argp;
1135:   struct DMCompositeLink *next;
1136:   DM_Composite           *com = (DM_Composite*)dm->data;

1140:   next = com->next;
1141:   /* loop over packed objects, handling one at at time */
1142:   va_start(Argp,dm);
1143:   while (next) {
1144:     DM *dmn;
1145:     dmn = va_arg(Argp, DM*);
1146:     if (dmn) *dmn = next->dm;
1147:     next = next->next;
1148:   }
1149:   va_end(Argp);
1150:   return(0);
1151: }

1155: /*@C
1156:     DMCompositeGetEntriesArray - Gets the DM for each entry in a DMComposite.

1158:     Not Collective

1160:     Input Parameter:
1161: .    dm - the packer object

1163:     Output Parameter:
1164: .    dms - array of sufficient length (see DMCompositeGetNumberDM()) to hold the individual DMs

1166:     Level: advanced

1168: .seealso DMDestroy(), DMCompositeAddDM(), DMCreateGlobalVector(), DMCompositeGetEntries()
1169:          DMCompositeGather(), DMCompositeCreate(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess(),
1170:          DMCompositeRestoreLocalVectors(), DMCompositeGetLocalVectors(),  DMCompositeScatter(),
1171:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors()

1173: @*/
1174: PetscErrorCode DMCompositeGetEntriesArray(DM dm,DM dms[])
1175: {
1176:   struct DMCompositeLink *next;
1177:   DM_Composite           *com = (DM_Composite*)dm->data;
1178:   PetscInt               i;

1182:   /* loop over packed objects, handling one at at time */
1183:   for (next=com->next,i=0; next; next=next->next,i++) dms[i] = next->dm;
1184:   return(0);
1185: }

1189: PetscErrorCode  DMRefine_Composite(DM dmi,MPI_Comm comm,DM *fine)
1190: {
1191:   PetscErrorCode         ierr;
1192:   struct DMCompositeLink *next;
1193:   DM_Composite           *com = (DM_Composite*)dmi->data;
1194:   DM                     dm;

1198:   if (comm == MPI_COMM_NULL) {
1199:     PetscObjectGetComm((PetscObject)dmi,&comm);
1200:   }
1201:   DMSetUp(dmi);
1202:   next = com->next;
1203:   DMCompositeCreate(comm,fine);

1205:   /* loop over packed objects, handling one at at time */
1206:   while (next) {
1207:     DMRefine(next->dm,comm,&dm);
1208:     DMCompositeAddDM(*fine,dm);
1209:     PetscObjectDereference((PetscObject)dm);
1210:     next = next->next;
1211:   }
1212:   return(0);
1213: }

1217: PetscErrorCode  DMCoarsen_Composite(DM dmi,MPI_Comm comm,DM *fine)
1218: {
1219:   PetscErrorCode         ierr;
1220:   struct DMCompositeLink *next;
1221:   DM_Composite           *com = (DM_Composite*)dmi->data;
1222:   DM                     dm;

1226:   DMSetUp(dmi);
1227:   if (comm == MPI_COMM_NULL) {
1228:     PetscObjectGetComm((PetscObject)dmi,&comm);
1229:   }
1230:   next = com->next;
1231:   DMCompositeCreate(comm,fine);

1233:   /* loop over packed objects, handling one at at time */
1234:   while (next) {
1235:     DMCoarsen(next->dm,comm,&dm);
1236:     DMCompositeAddDM(*fine,dm);
1237:     PetscObjectDereference((PetscObject)dm);
1238:     next = next->next;
1239:   }
1240:   return(0);
1241: }

1245: PetscErrorCode  DMCreateInterpolation_Composite(DM coarse,DM fine,Mat *A,Vec *v)
1246: {
1247:   PetscErrorCode         ierr;
1248:   PetscInt               m,n,M,N,nDM,i;
1249:   struct DMCompositeLink *nextc;
1250:   struct DMCompositeLink *nextf;
1251:   Vec                    gcoarse,gfine,*vecs;
1252:   DM_Composite           *comcoarse = (DM_Composite*)coarse->data;
1253:   DM_Composite           *comfine   = (DM_Composite*)fine->data;
1254:   Mat                    *mats;

1259:   DMSetUp(coarse);
1260:   DMSetUp(fine);
1261:   /* use global vectors only for determining matrix layout */
1262:   DMGetGlobalVector(coarse,&gcoarse);
1263:   DMGetGlobalVector(fine,&gfine);
1264:   VecGetLocalSize(gcoarse,&n);
1265:   VecGetLocalSize(gfine,&m);
1266:   VecGetSize(gcoarse,&N);
1267:   VecGetSize(gfine,&M);
1268:   DMRestoreGlobalVector(coarse,&gcoarse);
1269:   DMRestoreGlobalVector(fine,&gfine);

1271:   nDM = comfine->nDM;
1272:   if (nDM != comcoarse->nDM) SETERRQ2(PetscObjectComm((PetscObject)fine),PETSC_ERR_ARG_INCOMP,"Fine DMComposite has %D entries, but coarse has %D",nDM,comcoarse->nDM);
1273:   PetscCalloc1(nDM*nDM,&mats);
1274:   if (v) {
1275:     PetscCalloc1(nDM,&vecs);
1276:   }

1278:   /* loop over packed objects, handling one at at time */
1279:   for (nextc=comcoarse->next,nextf=comfine->next,i=0; nextc; nextc=nextc->next,nextf=nextf->next,i++) {
1280:     if (!v) {
1281:       DMCreateInterpolation(nextc->dm,nextf->dm,&mats[i*nDM+i],NULL);
1282:     } else {
1283:       DMCreateInterpolation(nextc->dm,nextf->dm,&mats[i*nDM+i],&vecs[i]);
1284:     }
1285:   }
1286:   MatCreateNest(PetscObjectComm((PetscObject)fine),nDM,NULL,nDM,NULL,mats,A);
1287:   if (v) {
1288:     VecCreateNest(PetscObjectComm((PetscObject)fine),nDM,NULL,vecs,v);
1289:   }
1290:   for (i=0; i<nDM*nDM; i++) {MatDestroy(&mats[i]);}
1291:   PetscFree(mats);
1292:   if (v) {
1293:     for (i=0; i<nDM; i++) {VecDestroy(&vecs[i]);}
1294:     PetscFree(vecs);
1295:   }
1296:   return(0);
1297: }

1301: static PetscErrorCode DMGetLocalToGlobalMapping_Composite(DM dm)
1302: {
1303:   DM_Composite           *com = (DM_Composite*)dm->data;
1304:   ISLocalToGlobalMapping *ltogs;
1305:   PetscInt               i;
1306:   PetscErrorCode         ierr;

1309:   /* Set the ISLocalToGlobalMapping on the new matrix */
1310:   DMCompositeGetISLocalToGlobalMappings(dm,&ltogs);
1311:   ISLocalToGlobalMappingConcatenate(PetscObjectComm((PetscObject)dm),com->nDM,ltogs,&dm->ltogmap);
1312:   for (i=0; i<com->nDM; i++) {ISLocalToGlobalMappingDestroy(&ltogs[i]);}
1313:   PetscFree(ltogs);
1314:   return(0);
1315: }


1320: PetscErrorCode  DMCreateColoring_Composite(DM dm,ISColoringType ctype,ISColoring *coloring)
1321: {
1322:   PetscErrorCode  ierr;
1323:   PetscInt        n,i,cnt;
1324:   ISColoringValue *colors;
1325:   PetscBool       dense  = PETSC_FALSE;
1326:   ISColoringValue maxcol = 0;
1327:   DM_Composite    *com   = (DM_Composite*)dm->data;

1331:   if (ctype == IS_COLORING_GHOSTED) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only global coloring supported");
1332:   else if (ctype == IS_COLORING_GLOBAL) {
1333:     n = com->n;
1334:   } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"Unknown ISColoringType");
1335:   PetscMalloc1(n,&colors); /* freed in ISColoringDestroy() */

1337:   PetscOptionsGetBool(NULL,"-dmcomposite_dense_jacobian",&dense,NULL);
1338:   if (dense) {
1339:     for (i=0; i<n; i++) {
1340:       colors[i] = (ISColoringValue)(com->rstart + i);
1341:     }
1342:     maxcol = com->N;
1343:   } else {
1344:     struct DMCompositeLink *next = com->next;
1345:     PetscMPIInt            rank;

1347:     MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);
1348:     cnt  = 0;
1349:     while (next) {
1350:       ISColoring lcoloring;

1352:       DMCreateColoring(next->dm,IS_COLORING_GLOBAL,&lcoloring);
1353:       for (i=0; i<lcoloring->N; i++) {
1354:         colors[cnt++] = maxcol + lcoloring->colors[i];
1355:       }
1356:       maxcol += lcoloring->n;
1357:       ISColoringDestroy(&lcoloring);
1358:       next    = next->next;
1359:     }
1360:   }
1361:   ISColoringCreate(PetscObjectComm((PetscObject)dm),maxcol,n,colors,PETSC_OWN_POINTER,coloring);
1362:   return(0);
1363: }

1367: PetscErrorCode  DMGlobalToLocalBegin_Composite(DM dm,Vec gvec,InsertMode mode,Vec lvec)
1368: {
1369:   PetscErrorCode         ierr;
1370:   struct DMCompositeLink *next;
1371:   PetscInt               cnt = 3;
1372:   PetscMPIInt            rank;
1373:   PetscScalar            *garray,*larray;
1374:   DM_Composite           *com = (DM_Composite*)dm->data;

1379:   next = com->next;
1380:   if (!com->setup) {
1381:     DMSetUp(dm);
1382:   }
1383:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);
1384:   VecGetArray(gvec,&garray);
1385:   VecGetArray(lvec,&larray);

1387:   /* loop over packed objects, handling one at at time */
1388:   while (next) {
1389:     Vec      local,global;
1390:     PetscInt N;

1392:     DMGetGlobalVector(next->dm,&global);
1393:     VecGetLocalSize(global,&N);
1394:     VecPlaceArray(global,garray);
1395:     DMGetLocalVector(next->dm,&local);
1396:     VecPlaceArray(local,larray);
1397:     DMGlobalToLocalBegin(next->dm,global,mode,local);
1398:     DMGlobalToLocalEnd(next->dm,global,mode,local);
1399:     VecResetArray(global);
1400:     VecResetArray(local);
1401:     DMRestoreGlobalVector(next->dm,&global);
1402:     DMRestoreGlobalVector(next->dm,&local);
1403:     cnt++;
1404:     larray += next->nlocal;
1405:     next    = next->next;
1406:   }

1408:   VecRestoreArray(gvec,NULL);
1409:   VecRestoreArray(lvec,NULL);
1410:   return(0);
1411: }

1415: PetscErrorCode  DMGlobalToLocalEnd_Composite(DM dm,Vec gvec,InsertMode mode,Vec lvec)
1416: {
1418:   return(0);
1419: }

1421: /*MC
1422:    DMCOMPOSITE = "composite" - A DM object that is used to manage data for a collection of DMs

1424:   Level: intermediate

1426: .seealso: DMType, DM, DMDACreate(), DMCreate(), DMSetType(), DMCompositeCreate()
1427: M*/


1432: PETSC_EXTERN PetscErrorCode DMCreate_Composite(DM p)
1433: {
1435:   DM_Composite   *com;

1438:   PetscNewLog(p,&com);
1439:   p->data   = com;
1440:   PetscObjectChangeTypeName((PetscObject)p,"DMComposite");
1441:   com->n    = 0;
1442:   com->next = NULL;
1443:   com->nDM  = 0;

1445:   p->ops->createglobalvector              = DMCreateGlobalVector_Composite;
1446:   p->ops->createlocalvector               = DMCreateLocalVector_Composite;
1447:   p->ops->getlocaltoglobalmapping         = DMGetLocalToGlobalMapping_Composite;
1448:   p->ops->createfieldis                   = DMCreateFieldIS_Composite;
1449:   p->ops->createfielddecomposition        = DMCreateFieldDecomposition_Composite;
1450:   p->ops->refine                          = DMRefine_Composite;
1451:   p->ops->coarsen                         = DMCoarsen_Composite;
1452:   p->ops->createinterpolation             = DMCreateInterpolation_Composite;
1453:   p->ops->creatematrix                    = DMCreateMatrix_Composite;
1454:   p->ops->getcoloring                     = DMCreateColoring_Composite;
1455:   p->ops->globaltolocalbegin              = DMGlobalToLocalBegin_Composite;
1456:   p->ops->globaltolocalend                = DMGlobalToLocalEnd_Composite;
1457:   p->ops->destroy                         = DMDestroy_Composite;
1458:   p->ops->view                            = DMView_Composite;
1459:   p->ops->setup                           = DMSetUp_Composite;
1460:   return(0);
1461: }

1465: /*@C
1466:     DMCompositeCreate - Creates a vector packer, used to generate "composite"
1467:       vectors made up of several subvectors.

1469:     Collective on MPI_Comm

1471:     Input Parameter:
1472: .   comm - the processors that will share the global vector

1474:     Output Parameters:
1475: .   packer - the packer object

1477:     Level: advanced

1479: .seealso DMDestroy(), DMCompositeAddDM(), DMCompositeScatter(), DMCOMPOSITE,DMCreate()
1480:          DMCompositeGather(), DMCreateGlobalVector(), DMCompositeGetISLocalToGlobalMappings(), DMCompositeGetAccess()
1481:          DMCompositeGetLocalVectors(), DMCompositeRestoreLocalVectors(), DMCompositeGetEntries()

1483: @*/
1484: PetscErrorCode  DMCompositeCreate(MPI_Comm comm,DM *packer)
1485: {

1490:   DMCreate(comm,packer);
1491:   DMSetType(*packer,DMCOMPOSITE);
1492:   return(0);
1493: }