Actual source code: dmshell.c

petsc-master 2020-01-26
Report Typos and Errors
  1:  #include <petscdmshell.h>
  2:  #include <petscmat.h>
  3:  #include <petsc/private/dmimpl.h>

  5: typedef struct  {
  6:   Vec        Xglobal;
  7:   Vec        Xlocal;
  8:   Mat        A;
  9:   VecScatter gtol;
 10:   VecScatter ltog;
 11:   VecScatter ltol;
 12:   void       *ctx;
 13: } DM_Shell;

 15: /*@
 16:    DMGlobalToLocalBeginDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to begin a global to local scatter
 17:    Collective

 19:    Input Arguments:
 20: +  dm - shell DM
 21: .  g - global vector
 22: .  mode - InsertMode
 23: -  l - local vector

 25:    Level: advanced

 27:    Note:  This is not normally called directly by user code, generally user code calls DMGlobalToLocalBegin() and DMGlobalToLocalEnd(). If the user provides their own custom routines to DMShellSetLocalToGlobal() then those routines might have reason to call this function.

 29: .seealso: DMGlobalToLocalEndDefaultShell()
 30: @*/
 31: PetscErrorCode DMGlobalToLocalBeginDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
 32: {
 34:   DM_Shell       *shell = (DM_Shell*)dm->data;

 37:   if (!shell->gtol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
 38:   VecScatterBegin(shell->gtol,g,l,mode,SCATTER_FORWARD);
 39:   return(0);
 40: }

 42: /*@
 43:    DMGlobalToLocalEndDefaultShell - Uses the GlobalToLocal VecScatter context set by the user to end a global to local scatter
 44:    Collective

 46:    Input Arguments:
 47: +  dm - shell DM
 48: .  g - global vector
 49: .  mode - InsertMode
 50: -  l - local vector

 52:    Level: advanced

 54: .seealso: DMGlobalToLocalBeginDefaultShell()
 55: @*/
 56: PetscErrorCode DMGlobalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
 57: {
 59:   DM_Shell       *shell = (DM_Shell*)dm->data;

 62:    if (!shell->gtol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
 63:   VecScatterEnd(shell->gtol,g,l,mode,SCATTER_FORWARD);
 64:   return(0);
 65: }

 67: /*@
 68:    DMLocalToGlobalBeginDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to begin a local to global scatter
 69:    Collective

 71:    Input Arguments:
 72: +  dm - shell DM
 73: .  l - local vector
 74: .  mode - InsertMode
 75: -  g - global vector

 77:    Level: advanced

 79:    Note:  This is not normally called directly by user code, generally user code calls DMLocalToGlobalBegin() and DMLocalToGlobalEnd(). If the user provides their own custom routines to DMShellSetLocalToGlobal() then those routines might have reason to call this function.

 81: .seealso: DMLocalToGlobalEndDefaultShell()
 82: @*/
 83: PetscErrorCode DMLocalToGlobalBeginDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
 84: {
 86:   DM_Shell       *shell = (DM_Shell*)dm->data;

 89:   if (!shell->ltog) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
 90:   VecScatterBegin(shell->ltog,l,g,mode,SCATTER_FORWARD);
 91:   return(0);
 92: }

 94: /*@
 95:    DMLocalToGlobalEndDefaultShell - Uses the LocalToGlobal VecScatter context set by the user to end a local to global scatter
 96:    Collective

 98:    Input Arguments:
 99: +  dm - shell DM
100: .  l - local vector
101: .  mode - InsertMode
102: -  g - global vector

104:    Level: advanced

106: .seealso: DMLocalToGlobalBeginDefaultShell()
107: @*/
108: PetscErrorCode DMLocalToGlobalEndDefaultShell(DM dm,Vec l,InsertMode mode,Vec g)
109: {
111:   DM_Shell       *shell = (DM_Shell*)dm->data;

114:    if (!shell->ltog) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToGlobalVecScatter()");
115:   VecScatterEnd(shell->ltog,l,g,mode,SCATTER_FORWARD);
116:   return(0);
117: }

119: /*@
120:    DMLocalToLocalBeginDefaultShell - Uses the LocalToLocal VecScatter context set by the user to begin a local to local scatter
121:    Collective

123:    Input Arguments:
124: +  dm - shell DM
125: .  g - the original local vector
126: -  mode - InsertMode

128:    Output Parameter:
129: .  l  - the local vector with correct ghost values

131:    Level: advanced

133:    Note:  This is not normally called directly by user code, generally user code calls DMLocalToLocalBegin() and DMLocalToLocalEnd(). If the user provides their own custom routines to DMShellSetLocalToLocal() then those routines might have reason to call this function.

135: .seealso: DMLocalToLocalEndDefaultShell()
136: @*/
137: PetscErrorCode DMLocalToLocalBeginDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
138: {
140:   DM_Shell       *shell = (DM_Shell*)dm->data;

143:   if (!shell->ltol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetLocalToLocalVecScatter()");
144:   VecScatterBegin(shell->ltol,g,l,mode,SCATTER_FORWARD);
145:   return(0);
146: }

148: /*@
149:    DMLocalToLocalEndDefaultShell - Uses the LocalToLocal VecScatter context set by the user to end a local to local scatter
150:    Collective

152:    Input Arguments:
153: +  dm - shell DM
154: .  g - the original local vector
155: -  mode - InsertMode

157:    Output Parameter:
158: .  l  - the local vector with correct ghost values

160:    Level: advanced

162: .seealso: DMLocalToLocalBeginDefaultShell()
163: @*/
164: PetscErrorCode DMLocalToLocalEndDefaultShell(DM dm,Vec g,InsertMode mode,Vec l)
165: {
167:   DM_Shell       *shell = (DM_Shell*)dm->data;

170:    if (!shell->ltol) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE, "Cannot be used without first setting the scatter context via DMShellSetGlobalToLocalVecScatter()");
171:   VecScatterEnd(shell->ltol,g,l,mode,SCATTER_FORWARD);
172:   return(0);
173: }

175: static PetscErrorCode DMCreateMatrix_Shell(DM dm,Mat *J)
176: {
178:   DM_Shell       *shell = (DM_Shell*)dm->data;
179:   Mat            A;

184:   if (!shell->A) {
185:     if (shell->Xglobal) {
186:       PetscInt m,M;
187:       PetscInfo(dm,"Naively creating matrix using global vector distribution without preallocation\n");
188:       VecGetSize(shell->Xglobal,&M);
189:       VecGetLocalSize(shell->Xglobal,&m);
190:       MatCreate(PetscObjectComm((PetscObject)dm),&shell->A);
191:       MatSetSizes(shell->A,m,m,M,M);
192:       MatSetType(shell->A,dm->mattype);
193:       MatSetUp(shell->A);
194:     } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
195:   }
196:   A = shell->A;
197:   /* the check below is tacky and incomplete */
198:   if (dm->mattype) {
199:     PetscBool flg,aij,seqaij,mpiaij;
200:     PetscObjectTypeCompare((PetscObject)A,dm->mattype,&flg);
201:     PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&seqaij);
202:     PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&mpiaij);
203:     PetscStrcmp(dm->mattype,MATAIJ,&aij);
204:     if (!flg) {
205:       if (!(aij && (seqaij || mpiaij))) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_NOTSAMETYPE,"Requested matrix of type %s, but only %s available",dm->mattype,((PetscObject)A)->type_name);
206:     }
207:   }
208:   if (((PetscObject)A)->refct < 2) { /* We have an exclusive reference so we can give it out */
209:     PetscObjectReference((PetscObject)A);
210:     MatZeroEntries(A);
211:     *J   = A;
212:   } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
213:     MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,J);
214:     MatZeroEntries(*J);
215:   }
216:   return(0);
217: }

219: PetscErrorCode DMCreateGlobalVector_Shell(DM dm,Vec *gvec)
220: {
222:   DM_Shell       *shell = (DM_Shell*)dm->data;
223:   Vec            X;

228:   *gvec = 0;
229:   X     = shell->Xglobal;
230:   if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()");
231:   if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
232:     PetscObjectReference((PetscObject)X);
233:     VecZeroEntries(X);
234:     *gvec = X;
235:   } else { /* Need to create a copy */
236:     VecDuplicate(X,gvec);
237:     VecZeroEntries(*gvec);
238:   }
239:   VecSetDM(*gvec,dm);
240:   return(0);
241: }

243: PetscErrorCode DMCreateLocalVector_Shell(DM dm,Vec *gvec)
244: {
246:   DM_Shell       *shell = (DM_Shell*)dm->data;
247:   Vec            X;

252:   *gvec = 0;
253:   X     = shell->Xlocal;
254:   if (!X) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetLocalVector() or DMShellSetCreateLocalVector()");
255:   if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */
256:     PetscObjectReference((PetscObject)X);
257:     VecZeroEntries(X);
258:     *gvec = X;
259:   } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
260:     VecDuplicate(X,gvec);
261:     VecZeroEntries(*gvec);
262:   }
263:   VecSetDM(*gvec,dm);
264:   return(0);
265: }

267: /*@
268:    DMShellSetContext - set some data to be usable by this DM

270:    Collective

272:    Input Arguments:
273: +  dm - shell DM
274: -  ctx - the context

276:    Level: advanced

278: .seealso: DMCreateMatrix(), DMShellGetContext()
279: @*/
280: PetscErrorCode DMShellSetContext(DM dm,void *ctx)
281: {
282:   DM_Shell       *shell = (DM_Shell*)dm->data;
284:   PetscBool      isshell;

288:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
289:   if (!isshell) return(0);
290:   shell->ctx = ctx;
291:   return(0);
292: }

294: /*@
295:    DMShellGetContext - Returns the user-provided context associated to the DM

297:    Collective

299:    Input Argument:
300: .  dm - shell DM

302:    Output Argument:
303: .  ctx - the context

305:    Level: advanced

307: .seealso: DMCreateMatrix(), DMShellSetContext()
308: @*/
309: PetscErrorCode DMShellGetContext(DM dm,void **ctx)
310: {
311:   DM_Shell       *shell = (DM_Shell*)dm->data;
313:   PetscBool      isshell;

317:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
318:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
319:   *ctx = shell->ctx;
320:   return(0);
321: }

323: /*@
324:    DMShellSetMatrix - sets a template matrix associated with the DMShell

326:    Collective

328:    Input Arguments:
329: +  dm - shell DM
330: -  J - template matrix

332:    Level: advanced

334: .seealso: DMCreateMatrix(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
335: @*/
336: PetscErrorCode DMShellSetMatrix(DM dm,Mat J)
337: {
338:   DM_Shell       *shell = (DM_Shell*)dm->data;
340:   PetscBool      isshell;

345:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
346:   if (!isshell) return(0);
347:   PetscObjectReference((PetscObject)J);
348:   MatDestroy(&shell->A);
349:   shell->A = J;
350:   return(0);
351: }

353: /*@C
354:    DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM

356:    Logically Collective on dm

358:    Input Arguments:
359: +  dm - the shell DM
360: -  func - the function to create a matrix

362:    Level: advanced

364: .seealso: DMCreateMatrix(), DMShellSetMatrix(), DMShellSetContext(), DMShellGetContext()
365: @*/
366: PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,Mat*))
367: {
370:   dm->ops->creatematrix = func;
371:   return(0);
372: }

374: /*@
375:    DMShellSetGlobalVector - sets a template global vector associated with the DMShell

377:    Logically Collective on dm

379:    Input Arguments:
380: +  dm - shell DM
381: -  X - template vector

383:    Level: advanced

385: .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector()
386: @*/
387: PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X)
388: {
389:   DM_Shell       *shell = (DM_Shell*)dm->data;
391:   PetscBool      isshell;
392:   DM             vdm;

397:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
398:   if (!isshell) return(0);
399:   VecGetDM(X,&vdm);
400:   /*
401:       if the vector proposed as the new base global vector for the DM is a DM vector associated
402:       with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
403:       we get a circular dependency that prevents the DM from being destroy when it should be.
404:       This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
405:       DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
406:       to set its input vector (which is associated with the DM) as the base global vector.
407:       Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
408:       for pointing out the problem.
409:    */
410:   if (vdm == dm) return(0);
411:   PetscObjectReference((PetscObject)X);
412:   VecDestroy(&shell->Xglobal);
413:   shell->Xglobal = X;
414:   return(0);
415: }

417: /*@C
418:    DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM

420:    Logically Collective

422:    Input Arguments:
423: +  dm - the shell DM
424: -  func - the creation routine

426:    Level: advanced

428: .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
429: @*/
430: PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
431: {
434:   dm->ops->createglobalvector = func;
435:   return(0);
436: }

438: /*@
439:    DMShellSetLocalVector - sets a template local vector associated with the DMShell

441:    Logically Collective on dm

443:    Input Arguments:
444: +  dm - shell DM
445: -  X - template vector

447:    Level: advanced

449: .seealso: DMCreateLocalVector(), DMShellSetMatrix(), DMShellSetCreateLocalVector()
450: @*/
451: PetscErrorCode DMShellSetLocalVector(DM dm,Vec X)
452: {
453:   DM_Shell       *shell = (DM_Shell*)dm->data;
455:   PetscBool      isshell;
456:   DM             vdm;

461:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
462:   if (!isshell) return(0);
463:   VecGetDM(X,&vdm);
464:   /*
465:       if the vector proposed as the new base global vector for the DM is a DM vector associated
466:       with the same DM then the current base global vector for the DM is ok and if we replace it with the new one
467:       we get a circular dependency that prevents the DM from being destroy when it should be.
468:       This occurs when SNESSet/GetNPC() is used with a SNES that does not have a user provided
469:       DM attached to it since the inner SNES (which shares the DM with the outer SNES) tries
470:       to set its input vector (which is associated with the DM) as the base global vector.
471:       Thanks to Juan P. Mendez Granado Re: [petsc-maint] Nonlinear conjugate gradien
472:       for pointing out the problem.
473:    */
474:   if (vdm == dm) return(0);
475:   PetscObjectReference((PetscObject)X);
476:   VecDestroy(&shell->Xlocal);
477:   shell->Xlocal = X;
478:   return(0);
479: }

481: /*@C
482:    DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the shell DM

484:    Logically Collective

486:    Input Arguments:
487: +  dm - the shell DM
488: -  func - the creation routine

490:    Level: advanced

492: .seealso: DMShellSetLocalVector(), DMShellSetCreateMatrix(), DMShellSetContext(), DMShellGetContext()
493: @*/
494: PetscErrorCode DMShellSetCreateLocalVector(DM dm,PetscErrorCode (*func)(DM,Vec*))
495: {
498:   dm->ops->createlocalvector = func;
499:   return(0);
500: }

502: /*@C
503:    DMShellSetGlobalToLocal - Sets the routines used to perform a global to local scatter

505:    Logically Collective on dm

507:    Input Arguments
508: +  dm - the shell DM
509: .  begin - the routine that begins the global to local scatter
510: -  end - the routine that ends the global to local scatter

512:    Notes:
513:     If these functions are not provided but DMShellSetGlobalToLocalVecScatter() is called then
514:    DMGlobalToLocalBeginDefaultShell()/DMGlobalToLocalEndDefaultShell() are used to to perform the transfers

516:    Level: advanced

518: .seealso: DMShellSetLocalToGlobal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
519: @*/
520: PetscErrorCode DMShellSetGlobalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
523:   dm->ops->globaltolocalbegin = begin;
524:   dm->ops->globaltolocalend = end;
525:   return(0);
526: }

528: /*@C
529:    DMShellSetLocalToGlobal - Sets the routines used to perform a local to global scatter

531:    Logically Collective on dm

533:    Input Arguments
534: +  dm - the shell DM
535: .  begin - the routine that begins the local to global scatter
536: -  end - the routine that ends the local to global scatter

538:    Notes:
539:     If these functions are not provided but DMShellSetLocalToGlobalVecScatter() is called then
540:    DMLocalToGlobalBeginDefaultShell()/DMLocalToGlobalEndDefaultShell() are used to to perform the transfers

542:    Level: advanced

544: .seealso: DMShellSetGlobalToLocal()
545: @*/
546: PetscErrorCode DMShellSetLocalToGlobal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
549:   dm->ops->localtoglobalbegin = begin;
550:   dm->ops->localtoglobalend = end;
551:   return(0);
552: }

554: /*@C
555:    DMShellSetLocalToLocal - Sets the routines used to perform a local to local scatter

557:    Logically Collective on dm

559:    Input Arguments
560: +  dm - the shell DM
561: .  begin - the routine that begins the local to local scatter
562: -  end - the routine that ends the local to local scatter

564:    Notes:
565:     If these functions are not provided but DMShellSetLocalToLocalVecScatter() is called then
566:    DMLocalToLocalBeginDefaultShell()/DMLocalToLocalEndDefaultShell() are used to to perform the transfers

568:    Level: advanced

570: .seealso: DMShellSetGlobalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
571: @*/
572: PetscErrorCode DMShellSetLocalToLocal(DM dm,PetscErrorCode (*begin)(DM,Vec,InsertMode,Vec),PetscErrorCode (*end)(DM,Vec,InsertMode,Vec)) {
575:   dm->ops->localtolocalbegin = begin;
576:   dm->ops->localtolocalend = end;
577:   return(0);
578: }

580: /*@
581:    DMShellSetGlobalToLocalVecScatter - Sets a VecScatter context for global to local communication

583:    Logically Collective on dm

585:    Input Arguments
586: +  dm - the shell DM
587: -  gtol - the global to local VecScatter context

589:    Level: advanced

591: .seealso: DMShellSetGlobalToLocal(), DMGlobalToLocalBeginDefaultShell(), DMGlobalToLocalEndDefaultShell()
592: @*/
593: PetscErrorCode DMShellSetGlobalToLocalVecScatter(DM dm, VecScatter gtol)
594: {
595:   DM_Shell       *shell = (DM_Shell*)dm->data;

601:   PetscObjectReference((PetscObject)gtol);
602:   VecScatterDestroy(&shell->gtol);
603:   shell->gtol = gtol;
604:   return(0);
605: }

607: /*@
608:    DMShellSetLocalToGlobalVecScatter - Sets a VecScatter context for local to global communication

610:    Logically Collective on dm

612:    Input Arguments
613: +  dm - the shell DM
614: -  ltog - the local to global VecScatter context

616:    Level: advanced

618: .seealso: DMShellSetLocalToGlobal(), DMLocalToGlobalBeginDefaultShell(), DMLocalToGlobalEndDefaultShell()
619: @*/
620: PetscErrorCode DMShellSetLocalToGlobalVecScatter(DM dm, VecScatter ltog)
621: {
622:   DM_Shell       *shell = (DM_Shell*)dm->data;

628:   PetscObjectReference((PetscObject)ltog);
629:   VecScatterDestroy(&shell->ltog);
630:   shell->ltog = ltog;
631:   return(0);
632: }

634: /*@
635:    DMShellSetLocalToLocalVecScatter - Sets a VecScatter context for local to local communication

637:    Logically Collective on dm

639:    Input Arguments
640: +  dm - the shell DM
641: -  ltol - the local to local VecScatter context

643:    Level: advanced

645: .seealso: DMShellSetLocalToLocal(), DMLocalToLocalBeginDefaultShell(), DMLocalToLocalEndDefaultShell()
646: @*/
647: PetscErrorCode DMShellSetLocalToLocalVecScatter(DM dm, VecScatter ltol)
648: {
649:   DM_Shell       *shell = (DM_Shell*)dm->data;

655:   PetscObjectReference((PetscObject)ltol);
656:   VecScatterDestroy(&shell->ltol);
657:   shell->ltol = ltol;
658:   return(0);
659: }

661: /*@C
662:    DMShellSetCoarsen - Set the routine used to coarsen the shell DM

664:    Logically Collective on dm

666:    Input Arguments
667: +  dm - the shell DM
668: -  coarsen - the routine that coarsens the DM

670:    Level: advanced

672: .seealso: DMShellSetRefine(), DMCoarsen(), DMShellGetCoarsen(), DMShellSetContext(), DMShellGetContext()
673: @*/
674: PetscErrorCode DMShellSetCoarsen(DM dm, PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*))
675: {
677:   PetscBool      isshell;

681:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
682:   if (!isshell) return(0);
683:   dm->ops->coarsen = coarsen;
684:   return(0);
685: }

687: /*@C
688:    DMShellGetCoarsen - Get the routine used to coarsen the shell DM

690:    Logically Collective on dm

692:    Input Argument:
693: .  dm - the shell DM

695:    Output Argument:
696: .  coarsen - the routine that coarsens the DM

698:    Level: advanced

700: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
701: @*/
702: PetscErrorCode DMShellGetCoarsen(DM dm, PetscErrorCode (**coarsen)(DM,MPI_Comm,DM*))
703: {
705:   PetscBool      isshell;

709:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
710:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
711:   *coarsen = dm->ops->coarsen;
712:   return(0);
713: }

715: /*@C
716:    DMShellSetRefine - Set the routine used to refine the shell DM

718:    Logically Collective on dm

720:    Input Arguments
721: +  dm - the shell DM
722: -  refine - the routine that refines the DM

724:    Level: advanced

726: .seealso: DMShellSetCoarsen(), DMRefine(), DMShellGetRefine(), DMShellSetContext(), DMShellGetContext()
727: @*/
728: PetscErrorCode DMShellSetRefine(DM dm, PetscErrorCode (*refine)(DM,MPI_Comm,DM*))
729: {
731:   PetscBool      isshell;

735:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
736:   if (!isshell) return(0);
737:   dm->ops->refine = refine;
738:   return(0);
739: }

741: /*@C
742:    DMShellGetRefine - Get the routine used to refine the shell DM

744:    Logically Collective on dm

746:    Input Argument:
747: .  dm - the shell DM

749:    Output Argument:
750: .  refine - the routine that refines the DM

752:    Level: advanced

754: .seealso: DMShellSetCoarsen(), DMCoarsen(), DMShellSetRefine(), DMRefine()
755: @*/
756: PetscErrorCode DMShellGetRefine(DM dm, PetscErrorCode (**refine)(DM,MPI_Comm,DM*))
757: {
759:   PetscBool      isshell;

763:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
764:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
765:   *refine = dm->ops->refine;
766:   return(0);
767: }

769: /*@C
770:    DMShellSetCreateInterpolation - Set the routine used to create the interpolation operator

772:    Logically Collective on dm

774:    Input Arguments
775: +  dm - the shell DM
776: -  interp - the routine to create the interpolation

778:    Level: advanced

780: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateInterpolation(), DMShellSetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
781: @*/
782: PetscErrorCode DMShellSetCreateInterpolation(DM dm, PetscErrorCode (*interp)(DM,DM,Mat*,Vec*))
783: {
785:   PetscBool      isshell;

789:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
790:   if (!isshell) return(0);
791:   dm->ops->createinterpolation = interp;
792:   return(0);
793: }

795: /*@C
796:    DMShellGetCreateInterpolation - Get the routine used to create the interpolation operator

798:    Logically Collective on dm

800:    Input Argument:
801: +  dm - the shell DM

803:    Output Argument:
804: -  interp - the routine to create the interpolation

806:    Level: advanced

808: .seealso: DMShellGetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
809: @*/
810: PetscErrorCode DMShellGetCreateInterpolation(DM dm, PetscErrorCode (**interp)(DM,DM,Mat*,Vec*))
811: {
813:   PetscBool      isshell;

817:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
818:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
819:   *interp = dm->ops->createinterpolation;
820:   return(0);
821: }

823: /*@C
824:    DMShellSetCreateRestriction - Set the routine used to create the restriction operator

826:    Logically Collective on dm

828:    Input Arguments
829: +  dm - the shell DM
830: -  striction- the routine to create the restriction

832:    Level: advanced

834: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellGetCreateRestriction(), DMShellSetContext(), DMShellGetContext()
835: @*/
836: PetscErrorCode DMShellSetCreateRestriction(DM dm, PetscErrorCode (*restriction)(DM,DM,Mat*))
837: {
839:   PetscBool      isshell;

843:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
844:   if (!isshell) return(0);
845:   dm->ops->createrestriction = restriction;
846:   return(0);
847: }

849: /*@C
850:    DMShellGetCreateRestriction - Get the routine used to create the restriction operator

852:    Logically Collective on dm

854:    Input Argument:
855: +  dm - the shell DM

857:    Output Argument:
858: -  restriction - the routine to create the restriction

860:    Level: advanced

862: .seealso: DMShellSetCreateInjection(), DMCreateInterpolation(), DMShellSetContext(), DMShellGetContext()
863: @*/
864: PetscErrorCode DMShellGetCreateRestriction(DM dm, PetscErrorCode (**restriction)(DM,DM,Mat*))
865: {
867:   PetscBool      isshell;

871:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
872:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
873:   *restriction = dm->ops->createrestriction;
874:   return(0);
875: }

877: /*@C
878:    DMShellSetCreateInjection - Set the routine used to create the injection operator

880:    Logically Collective on dm

882:    Input Arguments
883: +  dm - the shell DM
884: -  inject - the routine to create the injection

886:    Level: advanced

888: .seealso: DMShellSetCreateInterpolation(), DMCreateInjection(), DMShellGetCreateInjection(), DMShellSetContext(), DMShellGetContext()
889: @*/
890: PetscErrorCode DMShellSetCreateInjection(DM dm, PetscErrorCode (*inject)(DM,DM,Mat*))
891: {
893:   PetscBool      isshell;

897:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
898:   if (!isshell) return(0);
899:   dm->ops->createinjection = inject;
900:   return(0);
901: }

903: /*@C
904:    DMShellGetCreateInjection - Get the routine used to create the injection operator

906:    Logically Collective on dm

908:    Input Argument:
909: +  dm - the shell DM

911:    Output Argument:
912: -  inject - the routine to create the injection

914:    Level: advanced

916: .seealso: DMShellGetCreateInterpolation(), DMCreateInjection(), DMShellSetContext(), DMShellGetContext()
917: @*/
918: PetscErrorCode DMShellGetCreateInjection(DM dm, PetscErrorCode (**inject)(DM,DM,Mat*))
919: {
921:   PetscBool      isshell;

925:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
926:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
927:   *inject = dm->ops->createinjection;
928:   return(0);
929: }

931: /*@C
932:    DMShellSetCreateFieldDecomposition - Set the routine used to create a decomposition of fields for the shell DM

934:    Logically Collective on dm

936:    Input Arguments
937: +  dm - the shell DM
938: -  decomp - the routine to create the decomposition

940:    Level: advanced

942: .seealso: DMCreateFieldDecomposition(), DMShellSetContext(), DMShellGetContext()
943: @*/
944: PetscErrorCode DMShellSetCreateFieldDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,DM**))
945: {
947:   PetscBool      isshell;

951:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
952:   if (!isshell) return(0);
953:   dm->ops->createfielddecomposition = decomp;
954:   return(0);
955: }

957: /*@C
958:    DMShellSetCreateDomainDecomposition - Set the routine used to create a domain decomposition for the shell DM

960:    Logically Collective on dm

962:    Input Arguments
963: +  dm - the shell DM
964: -  decomp - the routine to create the decomposition

966:    Level: advanced

968: .seealso: DMCreateDomainDecomposition(), DMShellSetContext(), DMShellGetContext()
969: @*/
970: PetscErrorCode DMShellSetCreateDomainDecomposition(DM dm, PetscErrorCode (*decomp)(DM,PetscInt*,char***, IS**,IS**,DM**))
971: {
973:   PetscBool      isshell;

977:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
978:   if (!isshell) return(0);
979:   dm->ops->createdomaindecomposition = decomp;
980:   return(0);
981: }

983: /*@C
984:    DMShellSetCreateDomainDecompositionScatters - Set the routine used to create the scatter contexts for domain decomposition with a shell DM

986:    Logically Collective on dm

988:    Input Arguments
989: +  dm - the shell DM
990: -  scatter - the routine to create the scatters

992:    Level: advanced

994: .seealso: DMCreateDomainDecompositionScatters(), DMShellSetContext(), DMShellGetContext()
995: @*/
996: PetscErrorCode DMShellSetCreateDomainDecompositionScatters(DM dm, PetscErrorCode (*scatter)(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**))
997: {
999:   PetscBool      isshell;

1003:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1004:   if (!isshell) return(0);
1005:   dm->ops->createddscatters = scatter;
1006:   return(0);
1007: }

1009: /*@C
1010:    DMShellSetCreateSubDM - Set the routine used to create a sub DM from the shell DM

1012:    Logically Collective on dm

1014:    Input Arguments
1015: +  dm - the shell DM
1016: -  subdm - the routine to create the decomposition

1018:    Level: advanced

1020: .seealso: DMCreateSubDM(), DMShellGetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1021: @*/
1022: PetscErrorCode DMShellSetCreateSubDM(DM dm, PetscErrorCode (*subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1023: {
1025:   PetscBool      isshell;

1029:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1030:   if (!isshell) return(0);
1031:   dm->ops->createsubdm = subdm;
1032:   return(0);
1033: }

1035: /*@C
1036:    DMShellGetCreateSubDM - Get the routine used to create a sub DM from the shell DM

1038:    Logically Collective on dm

1040:    Input Argument:
1041: +  dm - the shell DM

1043:    Output Argument:
1044: -  subdm - the routine to create the decomposition

1046:    Level: advanced

1048: .seealso: DMCreateSubDM(), DMShellSetCreateSubDM(), DMShellSetContext(), DMShellGetContext()
1049: @*/
1050: PetscErrorCode DMShellGetCreateSubDM(DM dm, PetscErrorCode (**subdm)(DM,PetscInt,const PetscInt[],IS*,DM*))
1051: {
1053:   PetscBool      isshell;

1057:   PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);
1058:   if (!isshell) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only use with DMSHELL type DMs");
1059:   *subdm = dm->ops->createsubdm;
1060:   return(0);
1061: }

1063: static PetscErrorCode DMDestroy_Shell(DM dm)
1064: {
1066:   DM_Shell       *shell = (DM_Shell*)dm->data;

1069:   MatDestroy(&shell->A);
1070:   VecDestroy(&shell->Xglobal);
1071:   VecDestroy(&shell->Xlocal);
1072:   VecScatterDestroy(&shell->gtol);
1073:   VecScatterDestroy(&shell->ltog);
1074:   VecScatterDestroy(&shell->ltol);
1075:   /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
1076:   PetscFree(shell);
1077:   return(0);
1078: }

1080: static PetscErrorCode DMView_Shell(DM dm,PetscViewer v)
1081: {
1083:   DM_Shell       *shell = (DM_Shell*)dm->data;

1086:   VecView(shell->Xglobal,v);
1087:   return(0);
1088: }

1090: static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v)
1091: {
1093:   DM_Shell       *shell = (DM_Shell*)dm->data;

1096:   VecCreate(PetscObjectComm((PetscObject)dm),&shell->Xglobal);
1097:   VecLoad(shell->Xglobal,v);
1098:   return(0);
1099: }

1101: PetscErrorCode DMCreateSubDM_Shell(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1102: {

1106:   if (subdm) {DMShellCreate(PetscObjectComm((PetscObject) dm), subdm);}
1107:   DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1108:   return(0);
1109: }

1111: PETSC_EXTERN PetscErrorCode DMCreate_Shell(DM dm)
1112: {
1114:   DM_Shell       *shell;

1117:   PetscNewLog(dm,&shell);
1118:   dm->data = shell;

1120:   dm->ops->destroy            = DMDestroy_Shell;
1121:   dm->ops->createglobalvector = DMCreateGlobalVector_Shell;
1122:   dm->ops->createlocalvector  = DMCreateLocalVector_Shell;
1123:   dm->ops->creatematrix       = DMCreateMatrix_Shell;
1124:   dm->ops->view               = DMView_Shell;
1125:   dm->ops->load               = DMLoad_Shell;
1126:   dm->ops->globaltolocalbegin = DMGlobalToLocalBeginDefaultShell;
1127:   dm->ops->globaltolocalend   = DMGlobalToLocalEndDefaultShell;
1128:   dm->ops->localtoglobalbegin = DMLocalToGlobalBeginDefaultShell;
1129:   dm->ops->localtoglobalend   = DMLocalToGlobalEndDefaultShell;
1130:   dm->ops->localtolocalbegin  = DMLocalToLocalBeginDefaultShell;
1131:   dm->ops->localtolocalend    = DMLocalToLocalEndDefaultShell;
1132:   dm->ops->createsubdm        = DMCreateSubDM_Shell;
1133:   return(0);
1134: }

1136: /*@
1137:     DMShellCreate - Creates a shell DM object, used to manage user-defined problem data

1139:     Collective

1141:     Input Parameter:
1142: .   comm - the processors that will share the global vector

1144:     Output Parameters:
1145: .   shell - the shell DM

1147:     Level: advanced

1149: .seealso DMDestroy(), DMCreateGlobalVector(), DMCreateLocalVector(), DMShellSetContext(), DMShellGetContext()
1150: @*/
1151: PetscErrorCode  DMShellCreate(MPI_Comm comm,DM *dm)
1152: {

1157:   DMCreate(comm,dm);
1158:   DMSetType(*dm,DMSHELL);
1159:   DMSetUp(*dm);
1160:   return(0);
1161: }