Actual source code: dmts.c
 
   petsc-3.6.2 2015-10-02
   
  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,NULL,PETSC_FUNCTION);
 27:   PetscViewerBinaryRead(viewer,&kdm->ops->ifunctionview,1,NULL,PETSC_FUNCTION);
 28:   PetscViewerBinaryRead(viewer,&kdm->ops->ifunctionload,1,NULL,PETSC_FUNCTION);
 29:   if (kdm->ops->ifunctionload) {
 30:     (*kdm->ops->ifunctionload)(&kdm->ifunctionctx,viewer);
 31:   }
 32:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobian,1,NULL,PETSC_FUNCTION);
 33:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobianview,1,NULL,PETSC_FUNCTION);
 34:   PetscViewerBinaryRead(viewer,&kdm->ops->ijacobianload,1,NULL,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, DMTS_CLASSID, "DMTS", "DMTS", "DMTS", comm, DMTSDestroy, DMTSView);
115:   return(0);
116: }
120: /* Attaches the DMTS to the coarse level.
121:  * Under what conditions should we copy versus duplicate?
122:  */
123: static PetscErrorCode DMCoarsenHook_DMTS(DM dm,DM dmc,void *ctx)
124: {
128:   DMCopyDMTS(dm,dmc);
129:   return(0);
130: }
134: /* This could restrict auxiliary information to the coarse level.
135:  */
136: static PetscErrorCode DMRestrictHook_DMTS(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
137: {
140:   return(0);
141: }
145: static PetscErrorCode DMSubDomainHook_DMTS(DM dm,DM subdm,void *ctx)
146: {
150:   DMCopyDMTS(dm,subdm);
151:   return(0);
152: }
156: /* This could restrict auxiliary information to the coarse level.
157:  */
158: static PetscErrorCode DMSubDomainRestrictHook_DMTS(DM dm,VecScatter gscat,VecScatter lscat,DM subdm,void *ctx)
159: {
161:   return(0);
162: }
166: /*@C
167:    DMTSCopy - copies the information in a DMTS to another DMTS
169:    Not Collective
171:    Input Argument:
172: +  kdm - Original DMTS
173: -  nkdm - DMTS to receive the data, should have been created with DMTSCreate()
175:    Level: developer
177: .seealso: DMTSCreate(), DMTSDestroy()
178: @*/
179: PetscErrorCode DMTSCopy(DMTS kdm,DMTS nkdm)
180: {
186:   nkdm->ops->rhsfunction = kdm->ops->rhsfunction;
187:   nkdm->ops->rhsjacobian = kdm->ops->rhsjacobian;
188:   nkdm->ops->ifunction   = kdm->ops->ifunction;
189:   nkdm->ops->ijacobian   = kdm->ops->ijacobian;
190:   nkdm->ops->solution    = kdm->ops->solution;
191:   nkdm->ops->destroy     = kdm->ops->destroy;
192:   nkdm->ops->duplicate   = kdm->ops->duplicate;
194:   nkdm->rhsfunctionctx = kdm->rhsfunctionctx;
195:   nkdm->rhsjacobianctx = kdm->rhsjacobianctx;
196:   nkdm->ifunctionctx   = kdm->ifunctionctx;
197:   nkdm->ijacobianctx   = kdm->ijacobianctx;
198:   nkdm->solutionctx    = kdm->solutionctx;
200:   nkdm->data = kdm->data;
202:   /*
203:   nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
204:   nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
205:   nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
206:   */
208:   /* implementation specific copy hooks */
209:   if (kdm->ops->duplicate) {(*kdm->ops->duplicate)(kdm,nkdm);}
210:   return(0);
211: }
215: /*@C
216:    DMGetDMTS - get read-only private DMTS context from a DM
218:    Not Collective
220:    Input Argument:
221: .  dm - DM to be used with TS
223:    Output Argument:
224: .  tsdm - private DMTS context
226:    Level: developer
228:    Notes:
229:    Use DMGetDMTSWrite() if write access is needed. The DMTSSetXXX API should be used wherever possible.
231: .seealso: DMGetDMTSWrite()
232: @*/
233: PetscErrorCode DMGetDMTS(DM dm,DMTS *tsdm)
234: {
239:   *tsdm = (DMTS) dm->dmts;
240:   if (!*tsdm) {
241:     PetscInfo(dm,"Creating new DMTS\n");
242:     DMTSCreate(PetscObjectComm((PetscObject)dm),tsdm);
243:     dm->dmts = (PetscObject) *tsdm;
244:     DMCoarsenHookAdd(dm,DMCoarsenHook_DMTS,DMRestrictHook_DMTS,NULL);
245:     DMSubDomainHookAdd(dm,DMSubDomainHook_DMTS,DMSubDomainRestrictHook_DMTS,NULL);
246:   }
247:   return(0);
248: }
252: /*@C
253:    DMGetDMTSWrite - get write access to private DMTS context from a DM
255:    Not Collective
257:    Input Argument:
258: .  dm - DM to be used with TS
260:    Output Argument:
261: .  tsdm - private DMTS context
263:    Level: developer
265: .seealso: DMGetDMTS()
266: @*/
267: PetscErrorCode DMGetDMTSWrite(DM dm,DMTS *tsdm)
268: {
270:   DMTS           sdm;
274:   DMGetDMTS(dm,&sdm);
275:   if (!sdm->originaldm) sdm->originaldm = dm;
276:   if (sdm->originaldm != dm) {  /* Copy on write */
277:     DMTS oldsdm = sdm;
278:     PetscInfo(dm,"Copying DMTS due to write\n");
279:     DMTSCreate(PetscObjectComm((PetscObject)dm),&sdm);
280:     DMTSCopy(oldsdm,sdm);
281:     DMTSDestroy((DMTS*)&dm->dmts);
282:     dm->dmts = (PetscObject) sdm;
283:   }
284:   *tsdm = sdm;
285:   return(0);
286: }
290: /*@C
291:    DMCopyDMTS - copies a DM context to a new DM
293:    Logically Collective
295:    Input Arguments:
296: +  dmsrc - DM to obtain context from
297: -  dmdest - DM to add context to
299:    Level: developer
301:    Note:
302:    The context is copied by reference. This function does not ensure that a context exists.
304: .seealso: DMGetDMTS(), TSSetDM()
305: @*/
306: PetscErrorCode DMCopyDMTS(DM dmsrc,DM dmdest)
307: {
313:   DMTSDestroy((DMTS*)&dmdest->dmts);
314:   dmdest->dmts = dmsrc->dmts;
315:   PetscObjectReference(dmdest->dmts);
316:   DMCoarsenHookAdd(dmdest,DMCoarsenHook_DMTS,DMRestrictHook_DMTS,NULL);
317:   DMSubDomainHookAdd(dmdest,DMSubDomainHook_DMTS,DMSubDomainRestrictHook_DMTS,NULL);
318:   return(0);
319: }
323: /*@C
324:    DMTSSetIFunction - set TS implicit function evaluation function
326:    Not Collective
328:    Input Arguments:
329: +  dm - DM to be used with TS
330: .  func - function evaluation function, see TSSetIFunction() for calling sequence
331: -  ctx - context for residual evaluation
333:    Level: advanced
335:    Note:
336:    TSSetFunction() is normally used, but it calls this function internally because the user context is actually
337:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
338:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
340: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
341: @*/
342: PetscErrorCode DMTSSetIFunction(DM dm,TSIFunction func,void *ctx)
343: {
345:   DMTS           tsdm;
349:   DMGetDMTSWrite(dm,&tsdm);
350:   if (func) tsdm->ops->ifunction = func;
351:   if (ctx)  tsdm->ifunctionctx = ctx;
352:   return(0);
353: }
357: /*@C
358:    DMTSGetIFunction - get TS implicit residual evaluation function
360:    Not Collective
362:    Input Argument:
363: .  dm - DM to be used with TS
365:    Output Arguments:
366: +  func - function evaluation function, see TSSetIFunction() for calling sequence
367: -  ctx - context for residual evaluation
369:    Level: advanced
371:    Note:
372:    TSGetFunction() is normally used, but it calls this function internally because the user context is actually
373:    associated with the DM.
375: .seealso: DMTSSetContext(), DMTSSetFunction(), TSSetFunction()
376: @*/
377: PetscErrorCode DMTSGetIFunction(DM dm,TSIFunction *func,void **ctx)
378: {
380:   DMTS           tsdm;
384:   DMGetDMTS(dm,&tsdm);
385:   if (func) *func = tsdm->ops->ifunction;
386:   if (ctx)  *ctx = tsdm->ifunctionctx;
387:   return(0);
388: }
393: /*@C
394:    DMTSSetRHSFunction - set TS explicit residual evaluation function
396:    Not Collective
398:    Input Arguments:
399: +  dm - DM to be used with TS
400: .  func - RHS function evaluation function, see TSSetRHSFunction() for calling sequence
401: -  ctx - context for residual evaluation
403:    Level: advanced
405:    Note:
406:    TSSetRSHFunction() is normally used, but it calls this function internally because the user context is actually
407:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
408:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
410: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
411: @*/
412: PetscErrorCode DMTSSetRHSFunction(DM dm,TSRHSFunction func,void *ctx)
413: {
415:   DMTS           tsdm;
419:   DMGetDMTSWrite(dm,&tsdm);
420:   if (func) tsdm->ops->rhsfunction = func;
421:   if (ctx)  tsdm->rhsfunctionctx = ctx;
422:   return(0);
423: }
427: /*@C
428:    DMTSGetSolutionFunction - gets the TS solution evaluation function
430:    Not Collective
432:    Input Arguments:
433: .  dm - DM to be used with TS
435:    Output Parameters:
436: +  func - solution function evaluation function, see TSSetSolution() for calling sequence
437: -  ctx - context for solution evaluation
439:    Level: advanced
441: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
442: @*/
443: PetscErrorCode DMTSGetSolutionFunction(DM dm,TSSolutionFunction *func,void **ctx)
444: {
446:   DMTS           tsdm;
450:   DMGetDMTS(dm,&tsdm);
451:   if (func) *func = tsdm->ops->solution;
452:   if (ctx)  *ctx  = tsdm->solutionctx;
453:   return(0);
454: }
458: /*@C
459:    DMTSSetSolutionFunction - set TS solution evaluation function
461:    Not Collective
463:    Input Arguments:
464: +  dm - DM to be used with TS
465: .  func - solution function evaluation function, see TSSetSolution() for calling sequence
466: -  ctx - context for solution evaluation
468:    Level: advanced
470:    Note:
471:    TSSetSolutionFunction() is normally used, but it calls this function internally because the user context is actually
472:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
473:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
475: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
476: @*/
477: PetscErrorCode DMTSSetSolutionFunction(DM dm,TSSolutionFunction func,void *ctx)
478: {
480:   DMTS           tsdm;
484:   DMGetDMTSWrite(dm,&tsdm);
485:   if (func) tsdm->ops->solution = func;
486:   if (ctx)  tsdm->solutionctx   = ctx;
487:   return(0);
488: }
492: /*@C
493:    DMTSSetForcingFunction - set TS forcing function evaluation function
495:    Not Collective
497:    Input Arguments:
498: +  dm - DM to be used with TS
499: .  f - forcing function evaluation function; see TSForcingFunction
500: -  ctx - context for solution evaluation
502:    Level: advanced
504:    Note:
505:    TSSetForcingFunction() is normally used, but it calls this function internally because the user context is actually
506:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
507:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
509: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian(), TSSetForcingFunction(), DMTSGetForcingFunction()
510: @*/
511: PetscErrorCode DMTSSetForcingFunction(DM dm,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
512: {
514:   DMTS           tsdm;
518:   DMGetDMTSWrite(dm,&tsdm);
519:   if (f) tsdm->ops->forcing = f;
520:   if (ctx)  tsdm->forcingctx   = ctx;
521:   return(0);
522: }
527: /*@C
528:    DMTSGetForcingFunction - get TS forcing function evaluation function
530:    Not Collective
532:    Input Argument:
533: .   dm - DM to be used with TS
535:    Output Arguments:
536: +  f - forcing function evaluation function; see TSForcingFunction for details
537: -  ctx - context for solution evaluation
539:    Level: advanced
541:    Note:
542:    TSSetForcingFunction() is normally used, but it calls this function internally because the user context is actually
543:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
544:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
546: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian(), TSSetForcingFunction(), DMTSGetForcingFunction()
547: @*/
548: PetscErrorCode DMTSGetForcingFunction(DM dm,PetscErrorCode (**f)(TS,PetscReal,Vec,void*),void **ctx)
549: {
551:   DMTS           tsdm;
555:   DMGetDMTSWrite(dm,&tsdm);
556:   if (f) *f = tsdm->ops->forcing;
557:   if (ctx) *ctx = tsdm->forcingctx;
558:   return(0);
559: }
563: /*@C
564:    DMTSGetRHSFunction - get TS explicit residual evaluation function
566:    Not Collective
568:    Input Argument:
569: .  dm - DM to be used with TS
571:    Output Arguments:
572: +  func - residual evaluation function, see TSSetRHSFunction() for calling sequence
573: -  ctx - context for residual evaluation
575:    Level: advanced
577:    Note:
578:    TSGetFunction() is normally used, but it calls this function internally because the user context is actually
579:    associated with the DM.
581: .seealso: DMTSSetContext(), DMTSSetFunction(), TSSetFunction()
582: @*/
583: PetscErrorCode DMTSGetRHSFunction(DM dm,TSRHSFunction *func,void **ctx)
584: {
586:   DMTS           tsdm;
590:   DMGetDMTS(dm,&tsdm);
591:   if (func) *func = tsdm->ops->rhsfunction;
592:   if (ctx)  *ctx = tsdm->rhsfunctionctx;
593:   return(0);
594: }
598: /*@C
599:    DMTSSetIJacobian - set TS Jacobian evaluation function
601:    Not Collective
603:    Input Argument:
604: +  dm - DM to be used with TS
605: .  func - Jacobian evaluation function, see TSSetIJacobian() for calling sequence
606: -  ctx - context for residual evaluation
608:    Level: advanced
610:    Note:
611:    TSSetJacobian() is normally used, but it calls this function internally because the user context is actually
612:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
613:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
615: .seealso: DMTSSetContext(), TSSetFunction(), DMTSGetJacobian(), TSSetJacobian()
616: @*/
617: PetscErrorCode DMTSSetIJacobian(DM dm,TSIJacobian func,void *ctx)
618: {
620:   DMTS           sdm;
624:   DMGetDMTSWrite(dm,&sdm);
625:   if (func) sdm->ops->ijacobian = func;
626:   if (ctx)  sdm->ijacobianctx   = ctx;
627:   return(0);
628: }
632: /*@C
633:    DMTSGetIJacobian - get TS Jacobian evaluation function
635:    Not Collective
637:    Input Argument:
638: .  dm - DM to be used with TS
640:    Output Arguments:
641: +  func - Jacobian evaluation function, see TSSetIJacobian() for calling sequence
642: -  ctx - context for residual evaluation
644:    Level: advanced
646:    Note:
647:    TSGetJacobian() is normally used, but it calls this function internally because the user context is actually
648:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
649:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
651: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
652: @*/
653: PetscErrorCode DMTSGetIJacobian(DM dm,TSIJacobian *func,void **ctx)
654: {
656:   DMTS           tsdm;
660:   DMGetDMTS(dm,&tsdm);
661:   if (func) *func = tsdm->ops->ijacobian;
662:   if (ctx)  *ctx = tsdm->ijacobianctx;
663:   return(0);
664: }
669: /*@C
670:    DMTSSetRHSJacobian - set TS Jacobian evaluation function
672:    Not Collective
674:    Input Argument:
675: +  dm - DM to be used with TS
676: .  func - Jacobian evaluation function, see TSSetRHSJacobian() for calling sequence
677: -  ctx - context for residual evaluation
679:    Level: advanced
681:    Note:
682:    TSSetJacobian() is normally used, but it calls this function internally because the user context is actually
683:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
684:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
686: .seealso: DMTSSetContext(), TSSetFunction(), DMTSGetJacobian(), TSSetJacobian()
687: @*/
688: PetscErrorCode DMTSSetRHSJacobian(DM dm,TSRHSJacobian func,void *ctx)
689: {
691:   DMTS           tsdm;
695:   DMGetDMTSWrite(dm,&tsdm);
696:   if (func) tsdm->ops->rhsjacobian = func;
697:   if (ctx)  tsdm->rhsjacobianctx = ctx;
698:   return(0);
699: }
703: /*@C
704:    DMTSGetRHSJacobian - get TS Jacobian evaluation function
706:    Not Collective
708:    Input Argument:
709: .  dm - DM to be used with TS
711:    Output Arguments:
712: +  func - Jacobian evaluation function, see TSSetRHSJacobian() for calling sequence
713: -  ctx - context for residual evaluation
715:    Level: advanced
717:    Note:
718:    TSGetJacobian() is normally used, but it calls this function internally because the user context is actually
719:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
720:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
722: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
723: @*/
724: PetscErrorCode DMTSGetRHSJacobian(DM dm,TSRHSJacobian *func,void **ctx)
725: {
727:   DMTS           tsdm;
731:   DMGetDMTS(dm,&tsdm);
732:   if (func) *func = tsdm->ops->rhsjacobian;
733:   if (ctx)  *ctx = tsdm->rhsjacobianctx;
734:   return(0);
735: }
739: /*@C
740:    DMTSSetIFunctionSerialize - sets functions used to view and load a IFunction context
742:    Not Collective
744:    Input Arguments:
745: +  dm - DM to be used with TS
746: .  view - viewer function
747: -  load - loading function
749:    Level: advanced
751: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
752: @*/
753: PetscErrorCode DMTSSetIFunctionSerialize(DM dm,PetscErrorCode (*view)(void*,PetscViewer),PetscErrorCode (*load)(void**,PetscViewer))
754: {
756:   DMTS           tsdm;
760:   DMGetDMTSWrite(dm,&tsdm);
761:   tsdm->ops->ifunctionview = view;
762:   tsdm->ops->ifunctionload = load;
763:   return(0);
764: }
768: /*@C
769:    DMTSSetIJacobianSerialize - sets functions used to view and load a IJacobian context
771:    Not Collective
773:    Input Arguments:
774: +  dm - DM to be used with TS
775: .  view - viewer function
776: -  load - loading function
778:    Level: advanced
780: .seealso: DMTSSetContext(), TSSetFunction(), DMTSSetJacobian()
781: @*/
782: PetscErrorCode DMTSSetIJacobianSerialize(DM dm,PetscErrorCode (*view)(void*,PetscViewer),PetscErrorCode (*load)(void**,PetscViewer))
783: {
785:   DMTS           tsdm;
789:   DMGetDMTSWrite(dm,&tsdm);
790:   tsdm->ops->ijacobianview = view;
791:   tsdm->ops->ijacobianload = load;
792:   return(0);
793: }