Actual source code: dmts.c

petsc-3.5.4 2015-05-23
Report Typos and Errors
  1: #include <petsc-private/tsimpl.h>     /*I "petscts.h" I*/
  2: #include <petsc-private/dmimpl.h>

  6: static PetscErrorCode DMTSDestroy(DMTS *kdm)
  7: {

 11:   if (!*kdm) return(0);
 13:   if (--((PetscObject)(*kdm))->refct > 0) {*kdm = 0; return(0);}
 14:   if ((*kdm)->ops->destroy) {((*kdm)->ops->destroy)(*kdm);}
 15:   PetscHeaderDestroy(kdm);
 16:   return(0);
 17: }

 21: PetscErrorCode DMTSLoad(DMTS kdm,PetscViewer viewer)
 22: {

 26:   PetscViewerBinaryRead(viewer,&kdm->ops->ifunction,1,PETSC_FUNCTION);
 27:   PetscViewerBinaryRead(viewer,&kdm->ops->ifunctionview,1,PETSC_FUNCTION);
 28:   PetscViewerBinaryRead(viewer,&kdm->ops->ifunctionload,1,PETSC_FUNCTION);
 29:   if (kdm->ops->ifunctionload) {
 30:     (*kdm->ops->ifunctionload)(&kdm->ifunctionctx,viewer);
 31:   }
 32:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobian,1,PETSC_FUNCTION);
 33:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobianview,1,PETSC_FUNCTION);
 34:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobianload,1,PETSC_FUNCTION);
 35:   if (kdm->ops->ijacobianload) {
 36:     (*kdm->ops->ijacobianload)(&kdm->ijacobianctx,viewer);
 37:   }
 38:   return(0);
 39: }

 43: PetscErrorCode DMTSView(DMTS kdm,PetscViewer viewer)
 44: {
 46:   PetscBool      isascii,isbinary;

 49:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 50:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
 51:   if (isascii) {
 52: #if defined(PETSC_SERIALIZE_FUNCTIONS)
 53:     const char *fname;

 55:     PetscFPTFind(kdm->ops->ifunction,&fname);
 56:     if (fname) {
 57:       PetscViewerASCIIPrintf(viewer,"  IFunction used by TS: %s\n",fname);
 58:     }
 59:     PetscFPTFind(kdm->ops->ijacobian,&fname);
 60:     if (fname) {
 61:       PetscViewerASCIIPrintf(viewer,"  IJacobian function used by TS: %s\n",fname);
 62:     }
 63: #endif
 64:   } else if (isbinary) {
 65:     struct {
 66:       TSIFunction ifunction;
 67:     } funcstruct;
 68:     struct {
 69:       PetscErrorCode (*ifunctionview)(void*,PetscViewer);
 70:     } funcviewstruct;
 71:     struct {
 72:       PetscErrorCode (*ifunctionload)(void**,PetscViewer);
 73:     } funcloadstruct;
 74:     struct {
 75:       TSIJacobian ijacobian;
 76:     } jacstruct;
 77:     struct {
 78:       PetscErrorCode (*ijacobianview)(void*,PetscViewer);
 79:     } jacviewstruct;
 80:     struct {
 81:       PetscErrorCode (*ijacobianload)(void**,PetscViewer);
 82:     } jacloadstruct;

 84:     funcstruct.ifunction         = kdm->ops->ifunction;
 85:     funcviewstruct.ifunctionview = kdm->ops->ifunctionview;
 86:     funcloadstruct.ifunctionload = kdm->ops->ifunctionload;
 87:     PetscViewerBinaryWrite(viewer,&funcstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 88:     PetscViewerBinaryWrite(viewer,&funcviewstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 89:     PetscViewerBinaryWrite(viewer,&funcloadstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 90:     if (kdm->ops->ifunctionview) {
 91:       (*kdm->ops->ifunctionview)(kdm->ifunctionctx,viewer);
 92:     }
 93:     jacstruct.ijacobian = kdm->ops->ijacobian;
 94:     jacviewstruct.ijacobianview = kdm->ops->ijacobianview;
 95:     jacloadstruct.ijacobianload = kdm->ops->ijacobianload;
 96:     PetscViewerBinaryWrite(viewer,&jacstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 97:     PetscViewerBinaryWrite(viewer,&jacviewstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 98:     PetscViewerBinaryWrite(viewer,&jacloadstruct,1,PETSC_FUNCTION,PETSC_FALSE);
 99:     if (kdm->ops->ijacobianview) {
100:       (*kdm->ops->ijacobianview)(kdm->ijacobianctx,viewer);
101:     }
102:   }
103:   return(0);
104: }

108: static PetscErrorCode DMTSCreate(MPI_Comm comm,DMTS *kdm)
109: {

113:   TSInitializePackage();
114:   PetscHeaderCreate(*kdm, _p_DMTS, struct _DMTSOps, DMTS_CLASSID, "DMTS", "DMTS", "DMTS", comm, DMTSDestroy, DMTSView);
115:   PetscMemzero((*kdm)->ops, sizeof(struct _DMTSOps));
116:   return(0);
117: }

121: /* Attaches the DMTS to the coarse level.
122:  * Under what conditions should we copy versus duplicate?
123:  */
124: static PetscErrorCode DMCoarsenHook_DMTS(DM dm,DM dmc,void *ctx)
125: {

129:   DMCopyDMTS(dm,dmc);
130:   return(0);
131: }

135: /* This could restrict auxiliary information to the coarse level.
136:  */
137: static PetscErrorCode DMRestrictHook_DMTS(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
138: {

141:   return(0);
142: }

146: static PetscErrorCode DMSubDomainHook_DMTS(DM dm,DM subdm,void *ctx)
147: {

151:   DMCopyDMTS(dm,subdm);
152:   return(0);
153: }

157: /* This could restrict auxiliary information to the coarse level.
158:  */
159: static PetscErrorCode DMSubDomainRestrictHook_DMTS(DM dm,VecScatter gscat,VecScatter lscat,DM subdm,void *ctx)
160: {
162:   return(0);
163: }

167: /*@C
168:    DMTSCopy - copies the information in a DMTS to another DMTS

170:    Not Collective

172:    Input Argument:
173: +  kdm - Original DMTS
174: -  nkdm - DMTS to receive the data, should have been created with DMTSCreate()

176:    Level: developer

178: .seealso: DMTSCreate(), DMTSDestroy()
179: @*/
180: PetscErrorCode DMTSCopy(DMTS kdm,DMTS nkdm)
181: {

187:   nkdm->ops->rhsfunction = kdm->ops->rhsfunction;
188:   nkdm->ops->rhsjacobian = kdm->ops->rhsjacobian;
189:   nkdm->ops->ifunction   = kdm->ops->ifunction;
190:   nkdm->ops->ijacobian   = kdm->ops->ijacobian;
191:   nkdm->ops->solution    = kdm->ops->solution;
192:   nkdm->ops->destroy     = kdm->ops->destroy;
193:   nkdm->ops->duplicate   = kdm->ops->duplicate;

195:   nkdm->rhsfunctionctx = kdm->rhsfunctionctx;
196:   nkdm->rhsjacobianctx = kdm->rhsjacobianctx;
197:   nkdm->ifunctionctx   = kdm->ifunctionctx;
198:   nkdm->ijacobianctx   = kdm->ijacobianctx;
199:   nkdm->solutionctx    = kdm->solutionctx;

201:   nkdm->data = kdm->data;

203:   /*
204:   nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
205:   nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
206:   nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
207:   */

209:   /* implementation specific copy hooks */
210:   if (kdm->ops->duplicate) {(*kdm->ops->duplicate)(kdm,nkdm);}
211:   return(0);
212: }

216: /*@C
217:    DMGetDMTS - get read-only private DMTS context from a DM

219:    Not Collective

221:    Input Argument:
222: .  dm - DM to be used with TS

224:    Output Argument:
225: .  tsdm - private DMTS context

227:    Level: developer

229:    Notes:
230:    Use DMGetDMTSWrite() if write access is needed. The DMTSSetXXX API should be used wherever possible.

232: .seealso: DMGetDMTSWrite()
233: @*/
234: PetscErrorCode DMGetDMTS(DM dm,DMTS *tsdm)
235: {

240:   *tsdm = (DMTS) dm->dmts;
241:   if (!*tsdm) {
242:     PetscInfo(dm,"Creating new DMTS\n");
243:     DMTSCreate(PetscObjectComm((PetscObject)dm),tsdm);
244:     dm->dmts = (PetscObject) *tsdm;
245:     DMCoarsenHookAdd(dm,DMCoarsenHook_DMTS,DMRestrictHook_DMTS,NULL);
246:     DMSubDomainHookAdd(dm,DMSubDomainHook_DMTS,DMSubDomainRestrictHook_DMTS,NULL);
247:   }
248:   return(0);
249: }

253: /*@C
254:    DMGetDMTSWrite - get write access to private DMTS context from a DM

256:    Not Collective

258:    Input Argument:
259: .  dm - DM to be used with TS

261:    Output Argument:
262: .  tsdm - private DMTS context

264:    Level: developer

266: .seealso: DMGetDMTS()
267: @*/
268: PetscErrorCode DMGetDMTSWrite(DM dm,DMTS *tsdm)
269: {
271:   DMTS           sdm;

275:   DMGetDMTS(dm,&sdm);
276:   if (!sdm->originaldm) sdm->originaldm = dm;
277:   if (sdm->originaldm != dm) {  /* Copy on write */
278:     DMTS oldsdm = sdm;
279:     PetscInfo(dm,"Copying DMTS due to write\n");
280:     DMTSCreate(PetscObjectComm((PetscObject)dm),&sdm);
281:     DMTSCopy(oldsdm,sdm);
282:     DMTSDestroy((DMTS*)&dm->dmts);
283:     dm->dmts = (PetscObject) sdm;
284:   }
285:   *tsdm = sdm;
286:   return(0);
287: }

291: /*@C
292:    DMCopyDMTS - copies a DM context to a new DM

294:    Logically Collective

296:    Input Arguments:
297: +  dmsrc - DM to obtain context from
298: -  dmdest - DM to add context to

300:    Level: developer

302:    Note:
303:    The context is copied by reference. This function does not ensure that a context exists.

305: .seealso: DMGetDMTS(), TSSetDM()
306: @*/
307: PetscErrorCode DMCopyDMTS(DM dmsrc,DM dmdest)
308: {

314:   DMTSDestroy((DMTS*)&dmdest->dmts);
315:   dmdest->dmts = dmsrc->dmts;
316:   PetscObjectReference(dmdest->dmts);
317:   DMCoarsenHookAdd(dmdest,DMCoarsenHook_DMTS,DMRestrictHook_DMTS,NULL);
318:   DMSubDomainHookAdd(dmdest,DMSubDomainHook_DMTS,DMSubDomainRestrictHook_DMTS,NULL);
319:   return(0);
320: }

324: /*@C
325:    DMTSSetIFunction - set TS implicit function evaluation function

327:    Not Collective

329:    Input Arguments:
330: +  dm - DM to be used with TS
331: .  func - function evaluation function, see TSSetIFunction() for calling sequence
332: -  ctx - context for residual evaluation

334:    Level: advanced

336:    Note:
337:    TSSetFunction() is normally used, but it calls this function internally because the user context is actually
338:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
339:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

341: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
342: @*/
343: PetscErrorCode DMTSSetIFunction(DM dm,TSIFunction func,void *ctx)
344: {
346:   DMTS           tsdm;

350:   DMGetDMTSWrite(dm,&tsdm);
351:   if (func) tsdm->ops->ifunction = func;
352:   if (ctx)  tsdm->ifunctionctx = ctx;
353:   return(0);
354: }

358: /*@C
359:    DMTSGetIFunction - get TS implicit residual evaluation function

361:    Not Collective

363:    Input Argument:
364: .  dm - DM to be used with TS

366:    Output Arguments:
367: +  func - function evaluation function, see TSSetIFunction() for calling sequence
368: -  ctx - context for residual evaluation

370:    Level: advanced

372:    Note:
373:    TSGetFunction() is normally used, but it calls this function internally because the user context is actually
374:    associated with the DM.

376: .seealso: DMTSSetContext(), DMTSSetFunction(), TSSetFunction()
377: @*/
378: PetscErrorCode DMTSGetIFunction(DM dm,TSIFunction *func,void **ctx)
379: {
381:   DMTS           tsdm;

385:   DMGetDMTS(dm,&tsdm);
386:   if (func) *func = tsdm->ops->ifunction;
387:   if (ctx)  *ctx = tsdm->ifunctionctx;
388:   return(0);
389: }


394: /*@C
395:    DMTSSetRHSFunction - set TS explicit residual evaluation function

397:    Not Collective

399:    Input Arguments:
400: +  dm - DM to be used with TS
401: .  func - RHS function evaluation function, see TSSetRHSFunction() for calling sequence
402: -  ctx - context for residual evaluation

404:    Level: advanced

406:    Note:
407:    TSSetRSHFunction() is normally used, but it calls this function internally because the user context is actually
408:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
409:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

411: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
412: @*/
413: PetscErrorCode DMTSSetRHSFunction(DM dm,TSRHSFunction func,void *ctx)
414: {
416:   DMTS           tsdm;

420:   DMGetDMTSWrite(dm,&tsdm);
421:   if (func) tsdm->ops->rhsfunction = func;
422:   if (ctx)  tsdm->rhsfunctionctx = ctx;
423:   return(0);
424: }

428: /*@C
429:    DMTSGetSolutionFunction - gets the TS solution evaluation function

431:    Not Collective

433:    Input Arguments:
434: .  dm - DM to be used with TS

436:    Output Parameters:
437: +  func - solution function evaluation function, see TSSetSolution() for calling sequence
438: -  ctx - context for solution evaluation

440:    Level: advanced

442: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
443: @*/
444: PetscErrorCode DMTSGetSolutionFunction(DM dm,TSSolutionFunction *func,void **ctx)
445: {
447:   DMTS           tsdm;

451:   DMGetDMTS(dm,&tsdm);
452:   if (func) *func = tsdm->ops->solution;
453:   if (ctx)  *ctx  = tsdm->solutionctx;
454:   return(0);
455: }

459: /*@C
460:    DMTSSetSolutionFunction - set TS solution evaluation function

462:    Not Collective

464:    Input Arguments:
465: +  dm - DM to be used with TS
466: .  func - solution function evaluation function, see TSSetSolution() for calling sequence
467: -  ctx - context for solution evaluation

469:    Level: advanced

471:    Note:
472:    TSSetSolutionFunction() is normally used, but it calls this function internally because the user context is actually
473:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
474:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

476: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
477: @*/
478: PetscErrorCode DMTSSetSolutionFunction(DM dm,TSSolutionFunction func,void *ctx)
479: {
481:   DMTS           tsdm;

485:   DMGetDMTSWrite(dm,&tsdm);
486:   if (func) tsdm->ops->solution = func;
487:   if (ctx)  tsdm->solutionctx   = ctx;
488:   return(0);
489: }

493: /*@C
494:    DMTSSetForcingFunction - set TS forcing function evaluation function

496:    Not Collective

498:    Input Arguments:
499: +  dm - DM to be used with TS
500: .  f - forcing function evaluation function; see TSForcingFunction
501: -  ctx - context for solution evaluation

503:    Level: advanced

505:    Note:
506:    TSSetForcingFunction() is normally used, but it calls this function internally because the user context is actually
507:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
508:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

510: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian(), TSSetForcingFunction(), DMTSGetForcingFunction()
511: @*/
512: PetscErrorCode DMTSSetForcingFunction(DM dm,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
513: {
515:   DMTS           tsdm;

519:   DMGetDMTSWrite(dm,&tsdm);
520:   if (f) tsdm->ops->forcing = f;
521:   if (ctx)  tsdm->forcingctx   = ctx;
522:   return(0);
523: }


528: /*@C
529:    DMTSGetForcingFunction - get TS forcing function evaluation function

531:    Not Collective

533:    Input Argument:
534: .   dm - DM to be used with TS

536:    Output Arguments:
537: +  f - forcing function evaluation function; see TSForcingFunction for details
538: -  ctx - context for solution evaluation

540:    Level: advanced

542:    Note:
543:    TSSetForcingFunction() is normally used, but it calls this function internally because the user context is actually
544:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
545:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

547: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian(), TSSetForcingFunction(), DMTSGetForcingFunction()
548: @*/
549: PetscErrorCode DMTSGetForcingFunction(DM dm,PetscErrorCode (**f)(TS,PetscReal,Vec,void*),void **ctx)
550: {
552:   DMTS           tsdm;

556:   DMGetDMTSWrite(dm,&tsdm);
557:   if (f) *f = tsdm->ops->forcing;
558:   if (ctx) *ctx = tsdm->forcingctx;
559:   return(0);
560: }

564: /*@C
565:    DMTSGetRHSFunction - get TS explicit residual evaluation function

567:    Not Collective

569:    Input Argument:
570: .  dm - DM to be used with TS

572:    Output Arguments:
573: +  func - residual evaluation function, see TSSetRHSFunction() for calling sequence
574: -  ctx - context for residual evaluation

576:    Level: advanced

578:    Note:
579:    TSGetFunction() is normally used, but it calls this function internally because the user context is actually
580:    associated with the DM.

582: .seealso: DMTSSetContext(), DMTSSetFunction(), TSSetFunction()
583: @*/
584: PetscErrorCode DMTSGetRHSFunction(DM dm,TSRHSFunction *func,void **ctx)
585: {
587:   DMTS           tsdm;

591:   DMGetDMTS(dm,&tsdm);
592:   if (func) *func = tsdm->ops->rhsfunction;
593:   if (ctx)  *ctx = tsdm->rhsfunctionctx;
594:   return(0);
595: }

599: /*@C
600:    DMTSSetIJacobian - set TS Jacobian evaluation function

602:    Not Collective

604:    Input Argument:
605: +  dm - DM to be used with TS
606: .  func - Jacobian evaluation function, see TSSetIJacobian() for calling sequence
607: -  ctx - context for residual evaluation

609:    Level: advanced

611:    Note:
612:    TSSetJacobian() is normally used, but it calls this function internally because the user context is actually
613:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
614:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

616: .seealso: DMTSSetContext(), TSSetFunction(), DMTSGetJacobian(), TSSetJacobian()
617: @*/
618: PetscErrorCode DMTSSetIJacobian(DM dm,TSIJacobian func,void *ctx)
619: {
621:   DMTS           sdm;

625:   DMGetDMTSWrite(dm,&sdm);
626:   if (func) sdm->ops->ijacobian = func;
627:   if (ctx)  sdm->ijacobianctx   = ctx;
628:   return(0);
629: }

633: /*@C
634:    DMTSGetIJacobian - get TS Jacobian evaluation function

636:    Not Collective

638:    Input Argument:
639: .  dm - DM to be used with TS

641:    Output Arguments:
642: +  func - Jacobian evaluation function, see TSSetIJacobian() for calling sequence
643: -  ctx - context for residual evaluation

645:    Level: advanced

647:    Note:
648:    TSGetJacobian() is normally used, but it calls this function internally because the user context is actually
649:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
650:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

652: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
653: @*/
654: PetscErrorCode DMTSGetIJacobian(DM dm,TSIJacobian *func,void **ctx)
655: {
657:   DMTS           tsdm;

661:   DMGetDMTS(dm,&tsdm);
662:   if (func) *func = tsdm->ops->ijacobian;
663:   if (ctx)  *ctx = tsdm->ijacobianctx;
664:   return(0);
665: }


670: /*@C
671:    DMTSSetRHSJacobian - set TS Jacobian evaluation function

673:    Not Collective

675:    Input Argument:
676: +  dm - DM to be used with TS
677: .  func - Jacobian evaluation function, see TSSetRHSJacobian() for calling sequence
678: -  ctx - context for residual evaluation

680:    Level: advanced

682:    Note:
683:    TSSetJacobian() is normally used, but it calls this function internally because the user context is actually
684:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
685:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

687: .seealso: DMTSSetContext(), TSSetFunction(), DMTSGetJacobian(), TSSetJacobian()
688: @*/
689: PetscErrorCode DMTSSetRHSJacobian(DM dm,TSRHSJacobian func,void *ctx)
690: {
692:   DMTS           tsdm;

696:   DMGetDMTSWrite(dm,&tsdm);
697:   if (func) tsdm->ops->rhsjacobian = func;
698:   if (ctx)  tsdm->rhsjacobianctx = ctx;
699:   return(0);
700: }

704: /*@C
705:    DMTSGetRHSJacobian - get TS Jacobian evaluation function

707:    Not Collective

709:    Input Argument:
710: .  dm - DM to be used with TS

712:    Output Arguments:
713: +  func - Jacobian evaluation function, see TSSetRHSJacobian() for calling sequence
714: -  ctx - context for residual evaluation

716:    Level: advanced

718:    Note:
719:    TSGetJacobian() is normally used, but it calls this function internally because the user context is actually
720:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
721:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

723: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
724: @*/
725: PetscErrorCode DMTSGetRHSJacobian(DM dm,TSRHSJacobian *func,void **ctx)
726: {
728:   DMTS           tsdm;

732:   DMGetDMTS(dm,&tsdm);
733:   if (func) *func = tsdm->ops->rhsjacobian;
734:   if (ctx)  *ctx = tsdm->rhsjacobianctx;
735:   return(0);
736: }

740: /*@C
741:    DMTSSetIFunctionSerialize - sets functions used to view and load a IFunction context

743:    Not Collective

745:    Input Arguments:
746: +  dm - DM to be used with TS
747: .  view - viewer function
748: -  load - loading function

750:    Level: advanced

752: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
753: @*/
754: PetscErrorCode DMTSSetIFunctionSerialize(DM dm,PetscErrorCode (*view)(void*,PetscViewer),PetscErrorCode (*load)(void**,PetscViewer))
755: {
757:   DMTS           tsdm;

761:   DMGetDMTSWrite(dm,&tsdm);
762:   tsdm->ops->ifunctionview = view;
763:   tsdm->ops->ifunctionload = load;
764:   return(0);
765: }

769: /*@C
770:    DMTSSetIJacobianSerialize - sets functions used to view and load a IJacobian context

772:    Not Collective

774:    Input Arguments:
775: +  dm - DM to be used with TS
776: .  view - viewer function
777: -  load - loading function

779:    Level: advanced

781: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
782: @*/
783: PetscErrorCode DMTSSetIJacobianSerialize(DM dm,PetscErrorCode (*view)(void*,PetscViewer),PetscErrorCode (*load)(void**,PetscViewer))
784: {
786:   DMTS           tsdm;

790:   DMGetDMTSWrite(dm,&tsdm);
791:   tsdm->ops->ijacobianview = view;
792:   tsdm->ops->ijacobianload = load;
793:   return(0);
794: }