Actual source code: petscpvode.c

  1: /*$Id: petscpvode.c,v 1.68 2001/04/10 19:37:10 bsmith Exp bsmith $*/

 3:  #include petsc.h
  4: /*
  5:     Provides a PETSc interface to PVODE. Alan Hindmarsh's parallel ODE
  6:    solver.
  7: */

  9: #include "src/ts/impls/implicit/pvode/petscpvode.h"  /*I "petscts.h" I*/    

 11: /*
 12:       TSPrecond_PVode - function that we provide to PVODE to
 13:                         evaluate the preconditioner.

 15:     Contributed by: Liyang Xu

 17: */
 18: int TSPrecond_PVode(integer N,real tn,N_Vector y,N_Vector fy,boole jok,
 19:                     boole *jcurPtr,real _gamma,N_Vector ewt,real h,
 20:                     real uround,long int *nfePtr,void *P_data,
 21:                     N_Vector vtemp1,N_Vector vtemp2,N_Vector vtemp3)
 22: {
 23:   TS           ts = (TS) P_data;
 24:   TS_PVode     *cvode = (TS_PVode*)ts->data;
 25:   PC           pc = cvode->pc;
 26:   int          ierr;
 27:   Mat          Jac = ts->B;
 28:   Vec          tmpy = cvode->w1;
 29:   Scalar       one = 1.0,gm;
 30:   MatStructure str = DIFFERENT_NONZERO_PATTERN;
 31: 
 33:   /* This allows us to construct preconditioners in-place if we like */
 34:   MatSetUnfactored(Jac);

 36:   /*
 37:        jok - TRUE means reuse current Jacobian else recompute Jacobian
 38:   */
 39:   if (jok) {
 40:     ierr     = MatCopy(cvode->pmat,Jac,SAME_NONZERO_PATTERN);
 41:     str      = SAME_NONZERO_PATTERN;
 42:     *jcurPtr = FALSE;
 43:   } else {
 44:     /* make PETSc vector tmpy point to PVODE vector y */
 45:     VecPlaceArray(tmpy,&N_VIth(y,0));

 47:     /* compute the Jacobian */
 48:     TSComputeRHSJacobian(ts,ts->ptime,tmpy,&Jac,&Jac,&str);

 50:     /* copy the Jacobian matrix */
 51:     if (!cvode->pmat) {
 52:       MatDuplicate(Jac,MAT_COPY_VALUES,&cvode->pmat);
 53:       PetscLogObjectParent(ts,cvode->pmat);
 54:     }
 55:     MatCopy(Jac,cvode->pmat,SAME_NONZERO_PATTERN);

 57:     *jcurPtr = TRUE;
 58:   }

 60:   /* construct I-gamma*Jac  */
 61:   gm   = -_gamma;
 62:   MatScale(&gm,Jac);
 63:   MatShift(&one,Jac);
 64: 
 65:   PCSetOperators(pc,Jac,Jac,str);
 66:   return(0);
 67: }

 69: /*
 70:      TSPSolve_PVode -  routine that we provide to PVode that applies the preconditioner.
 71:       
 72:     Contributed by: Liyang Xu

 74: */
 75: int TSPSolve_PVode(integer N,real tn,N_Vector y,N_Vector fy,N_Vector vtemp,
 76:                    real _gamma,N_Vector ewt,real delta,long int *nfePtr,
 77:                    N_Vector r,int lr,void *P_data,N_Vector z)
 78: {
 79:   TS       ts = (TS) P_data;
 80:   TS_PVode *cvode = (TS_PVode*)ts->data;
 81:   PC       pc = cvode->pc;
 82:   Vec      rr = cvode->w1,xx = cvode->w2;
 83:   int      ierr;

 86:   /*
 87:       Make the PETSc work vectors rr and xx point to the arrays in the PVODE vectors 
 88:   */
 89:   VecPlaceArray(rr,&N_VIth(r,0));
 90:   VecPlaceArray(xx,&N_VIth(z,0));

 92:   /* 
 93:       Solve the Px=r and put the result in xx 
 94:   */
 95:   PCApply(pc,rr,xx);

 97:   return(0);
 98: }

100: /*
101:         TSFunction_PVode - routine that we provide to PVode that applies the right hand side.
102:       
103:     Contributed by: Liyang Xu
104: */
105: void TSFunction_PVode(int N,double t,N_Vector y,N_Vector ydot,void *ctx)
106: {
107:   TS        ts = (TS) ctx;
108:   TS_PVode *cvode = (TS_PVode*)ts->data;
109:   Vec       tmpx = cvode->w1,tmpy = cvode->w2;
110:   int       ierr;

113:   /*
114:       Make the PETSc work vectors tmpx and tmpy point to the arrays in the PVODE vectors 
115:   */
116:   VecPlaceArray(tmpx,&N_VIth(y,0));
117:   if (ierr) {
118:     (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
119:   }
120:   VecPlaceArray(tmpy,&N_VIth(ydot,0));
121:   if (ierr) {
122:     (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
123:   }

125:   /* now compute the right hand side function */
126:   TSComputeRHSFunction(ts,t,tmpx,tmpy);
127:   if (ierr) {
128:     (*PetscErrorPrintf)("TSFunction_PVode:Could not compute RHS function. Error code %d",ierr);
129:   }
130: }

132: /*
133:        TSStep_PVode_Nonlinear - Calls PVode to integrate the ODE.

135:     Contributed by: Liyang Xu
136: */
137: /* 
138:     TSStep_PVode_Nonlinear - 
139:   
140:    steps - number of time steps
141:    time - time that integrater is  terminated. 

143: */
144: int TSStep_PVode_Nonlinear(TS ts,int *steps,double *time)
145: {
146:   TS_PVode  *cvode = (TS_PVode*)ts->data;
147:   Vec       sol = ts->vec_sol;
148:   int       ierr,i,max_steps = ts->max_steps,flag;
149:   double    t,tout;

152:   /* initialize the number of steps */
153:   *steps = -ts->steps;
154:   ierr   = TSMonitor(ts,ts->steps,ts->ptime,sol);

156:   /* call CVSpgmr to use GMRES as the linear solver. */
157:   /* setup the ode integrator with the given preconditioner */
158:   CVSpgmr(cvode->mem,LEFT,cvode->gtype,cvode->restart,cvode->linear_tol,TSPrecond_PVode,TSPSolve_PVode,ts);

160:   tout = ts->max_time;
161:   for (i=0; i<max_steps; i++) {
162:     if (ts->ptime >= tout) break;
163:     VecGetArray(ts->vec_sol,&cvode->y->data);
164:     flag = CVode(cvode->mem,tout,cvode->y,&t,ONE_STEP);
165:     VecRestoreArray(ts->vec_sol,PETSC_NULL);
166:     if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE failed");

168:     if (t > tout && cvode->exact_final_time) {
169:       /* interpolate to final requested time */
170:       VecGetArray(ts->vec_sol,&cvode->y->data);
171:       flag = CVodeDky(cvode->mem,tout,0,cvode->y);
172:       VecRestoreArray(ts->vec_sol,PETSC_NULL);
173:       if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE interpolation to final time failed");
174:       t = tout;
175:     }

177:     ts->time_step = t - ts->ptime;
178:     ts->ptime     = t;

180:     /*
181:        copy the solution from cvode->y to cvode->update and sol 
182:     */
183:     VecPlaceArray(cvode->w1,&N_VIth(cvode->y,0));
184:     VecCopy(cvode->w1,cvode->update);
185:     VecCopy(cvode->update,sol);
186: 
187:     ts->steps++;
188:     TSMonitor(ts,ts->steps,t,sol);
189:     ts->nonlinear_its = cvode->iopt[NNI];
190:     ts->linear_its    = cvode->iopt[SPGMR_NLI];
191:   }

193:   *steps           += ts->steps;
194:   *time             = t;

196:   return(0);
197: }

199: /*

201:     Contributed by: Liyang Xu
202: */
203: int TSDestroy_PVode(TS ts)
204: {
205:   TS_PVode *cvode = (TS_PVode*)ts->data;
206:   int       ierr;

209:   if (cvode->pmat)   {MatDestroy(cvode->pmat);}
210:   if (cvode->pc)     {PCDestroy(cvode->pc);}
211:   if (cvode->update) {VecDestroy(cvode->update);}
212:   if (cvode->func)   {VecDestroy(cvode->func);}
213:   if (cvode->rhs)    {VecDestroy(cvode->rhs);}
214:   if (cvode->w1)     {VecDestroy(cvode->w1);}
215:   if (cvode->w2)     {VecDestroy(cvode->w2);}
216:   PetscFree(cvode);
217:   return(0);
218: }

220: /*

222:     Contributed by: Liyang Xu
223: */
224: int TSSetUp_PVode_Nonlinear(TS ts)
225: {
226:   TS_PVode    *cvode = (TS_PVode*)ts->data;
227:   int         ierr,M,locsize;
228:   machEnvType machEnv;

231:   PCSetFromOptions(cvode->pc);
232:   /* get the vector size */
233:   VecGetSize(ts->vec_sol,&M);
234:   VecGetLocalSize(ts->vec_sol,&locsize);

236:   /* allocate the memory for machEnv */
237:   /* machEnv = PVInitMPI(ts->comm,locsize,M);  */
238:   machEnv = PVecInitMPI(ts->comm, locsize, M, 0, 0);

240:   /* allocate the memory for N_Vec y */
241:   cvode->y         = N_VNew(M,machEnv);
242:   VecGetArray(ts->vec_sol,&cvode->y->data);
243:   VecRestoreArray(ts->vec_sol,PETSC_NULL);

245:   /* initializing vector update and func */
246:   VecDuplicate(ts->vec_sol,&cvode->update);
247:   VecDuplicate(ts->vec_sol,&cvode->func);
248:   PetscLogObjectParent(ts,cvode->update);
249:   PetscLogObjectParent(ts,cvode->func);

251:   /* 
252:       Create work vectors for the TSPSolve_PVode() routine. Note these are
253:     allocated with zero space arrays because the actual array space is provided 
254:     by PVode and set using VecPlaceArray().
255:   */
256:   VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w1);
257:   VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w2);
258:   PetscLogObjectParent(ts,cvode->w1);
259:   PetscLogObjectParent(ts,cvode->w2);

261:   PCSetVector(cvode->pc,ts->vec_sol);

263:   /* allocate memory for PVode */
264:   VecGetArray(ts->vec_sol,&cvode->y->data);
265:   cvode->mem = CVodeMalloc(M,TSFunction_PVode,ts->ptime,cvode->y,cvode->cvode_type,
266:                            NEWTON,SS,&cvode->reltol,&cvode->abstol,ts,NULL,FALSE,cvode->iopt,
267:                            cvode->ropt,machEnv);
268:   VecRestoreArray(ts->vec_sol,PETSC_NULL);
269:   return(0);
270: }

272: /*

274:     Contributed by: Liyang Xu
275: */
276: int TSSetFromOptions_PVode_Nonlinear(TS ts)
277: {
278:   TS_PVode   *cvode = (TS_PVode*)ts->data;
279:   int        ierr;
280:   char       method[128],*btype[] = {"bdf","adams"},*otype[] = {"modified","unmodified"};
281:   PetscTruth flag;

284:   PetscOptionsHead("PVODE ODE solver options");
285:     PetscOptionsEList("-ts_pvode_type","Scheme","TSPVodeSetType",btype,2,"bdf",method,127,&flag);
286:     if (flag) {
287:       PetscTruth isbdf,isadams;

289:       PetscStrcmp(method,btype[0],&isbdf);
290:       PetscStrcmp(method,btype[1],&isadams);
291:       if (isbdf) {
292:         TSPVodeSetType(ts,PVODE_BDF);
293:       } else if (isadams) {
294:         TSPVodeSetType(ts,PVODE_ADAMS);
295:       } else {
296:         SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode type.n");
297:       }
298:     }
299:     PetscOptionsEList("-ts_pvode_gramschmidt_type","Type of orthogonalization","TSPVodeSetGramSchmidtType",otype,2,"unmodified",method,127,&flag);
300:     if (flag) {
301:       PetscTruth ismodified,isunmodified;

303:       PetscStrcmp(method,otype[0],&ismodified);
304:       PetscStrcmp(method,otype[1],&isunmodified);
305:       if (ismodified) {
306:         TSPVodeSetGramSchmidtType(ts,PVODE_MODIFIED_GS);
307:       } else if (isunmodified) {
308:         TSPVodeSetGramSchmidtType(ts,PVODE_UNMODIFIED_GS);
309:       } else {
310:         SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode Gram-Schmidt orthogonalization type n");
311:       }
312:     }
313:     PetscOptionsDouble("-ts_pvode_atol","Absolute tolerance for convergence","TSPVodeSetTolerance",cvode->abstol,&cvode->abstol,PETSC_NULL);
314:     PetscOptionsDouble("-ts_pvode_rtol","Relative tolerance for convergence","TSPVodeSetTolerance",cvode->reltol,&cvode->reltol,PETSC_NULL);
315:     PetscOptionsDouble("-ts_pvode_linear_tolerance","Convergence tolerance for linear solve","TSPVodeSetLinearTolerance",cvode->linear_tol,&cvode->linear_tol,&flag);
316:     PetscOptionsInt("-ts_pvode_gmres_restart","Number of GMRES orthogonalization directions","TSPVodeSetGMRESRestart",cvode->restart,&cvode->restart,&flag);
317:     PetscOptionsName("-ts_pvode_not_exact_final_time","Allow PVODE to stop near the final time, not exactly on it","TSPVodeSetExactFinalTime",&cvode->exact_final_time);
318:   PetscOptionsTail();

320:   return(0);
321: }

323: /*

325:     Contributed by: Liyang Xu
326: */
327: int TSPrintHelp_PVode(TS ts,char *p)
328: {
329:   int      ierr;

332:   (*PetscHelpPrintf)(ts->comm," Options for TSPVODE integrater:n");
333:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_type <bdf,adams>: integration approachn",p);
334:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_atol aabs: absolute tolerance of ODE solutionn",p);
335:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_rtol rel: relative tolerance of ODE solutionn",p);
336:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_gramschmidt_type <unmodified,modified>n");
337:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_gmres_restart <restart_size> (also max. GMRES its)n");
338:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_linear_tolerance <tol>n");
339:   (*PetscHelpPrintf)(ts->comm," -ts_pvode_not_exact_final_timen");

341:   return(0);
342: }

344: /*

346:     Contributed by: Liyang Xu
347: */
348: int TSView_PVode(TS ts,PetscViewer viewer)
349: {
350:   TS_PVode   *cvode = (TS_PVode*)ts->data;
351:   int        ierr;
352:   char       *type;
353:   PetscTruth isascii,isstring;

356:   if (cvode->cvode_type == PVODE_ADAMS) {type = "Adams";}
357:   else                                  {type = "BDF: backward differentiation formula";}

359:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
360:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
361:   if (isascii) {
362:     PetscViewerASCIIPrintf(viewer,"PVode integrater does not use SNES!n");
363:     PetscViewerASCIIPrintf(viewer,"PVode integrater type %sn",type);
364:     PetscViewerASCIIPrintf(viewer,"PVode abs tol %g rel tol %gn",cvode->abstol,cvode->reltol);
365:     PetscViewerASCIIPrintf(viewer,"PVode linear solver tolerance factor %gn",cvode->linear_tol);
366:     PetscViewerASCIIPrintf(viewer,"PVode GMRES max iterations (same as restart in PVODE) %dn",cvode->restart);
367:     if (cvode->gtype == PVODE_MODIFIED_GS) {
368:       PetscViewerASCIIPrintf(viewer,"PVode using modified Gram-Schmidt for orthogonalization in GMRESn");
369:     } else {
370:       PetscViewerASCIIPrintf(viewer,"PVode using unmodified (classical) Gram-Schmidt for orthogonalization in GMRESn");
371:     }
372:   } else if (isstring) {
373:     PetscViewerStringSPrintf(viewer,"Pvode type %s",type);
374:   } else {
375:     SETERRQ1(1,"Viewer type %s not supported by TS PVode",((PetscObject)viewer)->type_name);
376:   }
377:   PetscViewerASCIIPushTab(viewer);
378:   PCView(cvode->pc,viewer);
379:   PetscViewerASCIIPopTab(viewer);

381:   return(0);
382: }


385: /* --------------------------------------------------------------------------*/
386: EXTERN_C_BEGIN
387: int TSPVodeSetType_PVode(TS ts,TSPVodeType type)
388: {
389:   TS_PVode *cvode = (TS_PVode*)ts->data;
390: 
392:   cvode->cvode_type = type;
393:   return(0);
394: }
395: EXTERN_C_END

397: EXTERN_C_BEGIN
398: int TSPVodeSetGMRESRestart_PVode(TS ts,int restart)
399: {
400:   TS_PVode *cvode = (TS_PVode*)ts->data;
401: 
403:   cvode->restart = restart;
404:   return(0);
405: }
406: EXTERN_C_END

408: EXTERN_C_BEGIN
409: int TSPVodeSetLinearTolerance_PVode(TS ts,double tol)
410: {
411:   TS_PVode *cvode = (TS_PVode*)ts->data;
412: 
414:   cvode->linear_tol = tol;
415:   return(0);
416: }
417: EXTERN_C_END

419: EXTERN_C_BEGIN
420: int TSPVodeSetGramSchmidtType_PVode(TS ts,TSPVodeGramSchmidtType type)
421: {
422:   TS_PVode *cvode = (TS_PVode*)ts->data;
423: 
425:   cvode->gtype = type;

427:   return(0);
428: }
429: EXTERN_C_END

431: EXTERN_C_BEGIN
432: int TSPVodeSetTolerance_PVode(TS ts,double aabs,double rel)
433: {
434:   TS_PVode *cvode = (TS_PVode*)ts->data;
435: 
437:   if (aabs != PETSC_DECIDE) cvode->abstol = aabs;
438:   if (rel != PETSC_DECIDE)  cvode->reltol = rel;
439:   return(0);
440: }
441: EXTERN_C_END

443: EXTERN_C_BEGIN
444: int TSPVodeGetPC_PVode(TS ts,PC *pc)
445: {
446:   TS_PVode *cvode = (TS_PVode*)ts->data;

449:   *pc = cvode->pc;

451:   return(0);
452: }
453: EXTERN_C_END

455: EXTERN_C_BEGIN
456: int TSPVodeGetIterations_PVode(TS ts,int *nonlin,int *lin)
457: {
458:   TS_PVode *cvode = (TS_PVode*)ts->data;
459: 
461:   if (nonlin) *nonlin = cvode->iopt[NNI];
462:   if (lin)    *lin    = cvode->iopt[SPGMR_NLI];
463:   return(0);
464: }
465: EXTERN_C_END
466:   
467: EXTERN_C_BEGIN
468: int TSPVodeSetExactFinalTime_PVode(TS ts,PetscTruth s)
469: {
470:   TS_PVode *cvode = (TS_PVode*)ts->data;
471: 
473:   cvode->exact_final_time = s;
474:   return(0);
475: }
476: EXTERN_C_END
477: /* -------------------------------------------------------------------------------------------*/

479: /*@C
480:    TSPVodeGetIterations - Gets the number of nonlinear and linear iterations used so far by PVode.

482:    Not Collective

484:    Input parameters:
485: .    ts     - the time-step context

487:    Output Parameters:
488: +   nonlin - number of nonlinear iterations
489: -   lin    - number of linear iterations

491:    Level: advanced

493: .keywords: non-linear iterations, linear iterations

495: .seealso: TSPVodeSetType(), TSPVodeSetGMRESRestart(),
496:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
497:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
498:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
499:           TSPVodeSetExactFinalTime()

501: @*/
502: int TSPVodeGetIterations(TS ts,int *nonlin,int *lin)
503: {
504:   int ierr,(*f)(TS,int*,int*);
505: 
507:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetIterations_C",(void (**)())&f);
508:   if (f) {
509:     (*f)(ts,nonlin,lin);
510:   }
511:   return(0);
512: }

514: /*@
515:    TSPVodeSetType - Sets the method that PVode will use for integration.

517:    Collective on TS

519:    Input parameters:
520: +    ts     - the time-step context
521: -    type   - one of  PVODE_ADAMS or PVODE_BDF

523:     Contributed by: Liyang Xu

525:    Level: intermediate

527: .keywords: Adams, backward differentiation formula

529: .seealso: TSPVodeGetIterations(),  TSPVodeSetGMRESRestart(),
530:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
531:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
532:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
533:           TSPVodeSetExactFinalTime()
534: @*/
535: int TSPVodeSetType(TS ts,TSPVodeType type)
536: {
537:   int ierr,(*f)(TS,TSPVodeType);
538: 
540:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetType_C",(void (**)())&f);
541:   if (f) {
542:     (*f)(ts,type);
543:   }
544:   return(0);
545: }

547: /*@
548:    TSPVodeSetGMRESRestart - Sets the dimension of the Krylov space used by 
549:        GMRES in the linear solver in PVODE. PVODE DOES NOT use restarted GMRES so
550:        this is ALSO the maximum number of GMRES steps that will be used.

552:    Collective on TS

554:    Input parameters:
555: +    ts      - the time-step context
556: -    restart - number of direction vectors (the restart size).

558:    Level: advanced

560: .keywords: GMRES, restart

562: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), 
563:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
564:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
565:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
566:           TSPVodeSetExactFinalTime()

568: @*/
569: int TSPVodeSetGMRESRestart(TS ts,int restart)
570: {
571:   int ierr,(*f)(TS,int);

574:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGMRESRestart_C",(void (**)())&f);
575:   if (f) {
576:     (*f)(ts,restart);
577:   }

579:   return(0);
580: }

582: /*@
583:    TSPVodeSetLinearTolerance - Sets the tolerance used to solve the linear
584:        system by PVODE.

586:    Collective on TS

588:    Input parameters:
589: +    ts     - the time-step context
590: -    tol    - the factor by which the tolerance on the nonlinear solver is
591:              multiplied to get the tolerance on the linear solver, .05 by default.

593:    Level: advanced

595: .keywords: GMRES, linear convergence tolerance, PVODE

597: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
598:           TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
599:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
600:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
601:           TSPVodeSetExactFinalTime()

603: @*/
604: int TSPVodeSetLinearTolerance(TS ts,double tol)
605: {
606:   int ierr,(*f)(TS,double);
607: 
609:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetLinearTolerance_C",(void (**)())&f);
610:   if (f) {
611:     (*f)(ts,tol);
612:   }
613:   return(0);
614: }

616: /*@
617:    TSPVodeSetGramSchmidtType - Sets type of orthogonalization used
618:         in GMRES method by PVODE linear solver.

620:    Collective on TS

622:    Input parameters:
623: +    ts  - the time-step context
624: -    type - either PVODE_MODIFIED_GS or PVODE_CLASSICAL_GS

626:    Level: advanced

628: .keywords: PVode, orthogonalization

630: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
631:           TSPVodeSetLinearTolerance(),  TSPVodeSetTolerance(),
632:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
633:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
634:           TSPVodeSetExactFinalTime()

636: @*/
637: int TSPVodeSetGramSchmidtType(TS ts,TSPVodeGramSchmidtType type)
638: {
639:   int ierr,(*f)(TS,TSPVodeGramSchmidtType);
640: 
642:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",(void (**)())&f);
643:   if (f) {
644:     (*f)(ts,type);
645:   }
646:   return(0);
647: }

649: /*@
650:    TSPVodeSetTolerance - Sets the absolute and relative tolerance used by 
651:                          PVode for error control.

653:    Collective on TS

655:    Input parameters:
656: +    ts  - the time-step context
657: .    aabs - the absolute tolerance  
658: -    rel - the relative tolerance

660:     Contributed by: Liyang Xu

662:    Level: intermediate

664: .keywords: PVode, tolerance

666: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
667:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), 
668:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
669:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
670:           TSPVodeSetExactFinalTime()

672: @*/
673: int TSPVodeSetTolerance(TS ts,double aabs,double rel)
674: {
675:   int ierr,(*f)(TS,double,double);
676: 
678:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetTolerance_C",(void (**)())&f);
679:   if (f) {
680:     (*f)(ts,aabs,rel);
681:   }
682:   return(0);
683: }

685: /*@
686:    TSPVodeGetPC - Extract the PC context from a time-step context for PVode.

688:    Input Parameter:
689: .    ts - the time-step context

691:    Output Parameter:
692: .    pc - the preconditioner context

694:    Level: advanced

696:     Contributed by: Liyang Xu

698: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
699:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
700:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
701:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance()
702: @*/
703: int TSPVodeGetPC(TS ts,PC *pc)
704: {
705:   int ierr,(*f)(TS,PC *);

708:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetPC_C",(void (**)())&f);
709:   if (f) {
710:     (*f)(ts,pc);
711:   } else {
712:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"TS must be of PVode type to extract the PC");
713:   }

715:   return(0);
716: }

718: /*@
719:    TSPVodeSetExactFinalTime - Determines if PVode interpolates solution to the 
720:       exact final time requested by the user or just returns it at the final time
721:       it computed. (Defaults to true).

723:    Input Parameter:
724: +   ts - the time-step context
725: -   ft - PETSC_TRUE if interpolates, else PETSC_FALSE

727:    Level: beginner

729: .seealso:TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
730:           TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
731:           TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
732:           TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC() 
733: @*/
734: int TSPVodeSetExactFinalTime(TS ts,PetscTruth ft)
735: {
736:   int ierr,(*f)(TS,PetscTruth);

739:   PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetExactFinalTime_C",(void (**)())&f);
740:   if (f) {
741:     (*f)(ts,ft);
742:   }

744:   return(0);
745: }

747: /* -------------------------------------------------------------------------------------------*/

749: /*

751:     Contributed by: Liyang Xu
752: */
753: EXTERN_C_BEGIN
754: int TSCreate_PVode(TS ts)
755: {
756:   TS_PVode *cvode;
757:   int      ierr;

760:   ts->destroy         = TSDestroy_PVode;
761:   ts->view            = TSView_PVode;

763:   if (ts->problem_type != TS_NONLINEAR) {
764:     SETERRQ(PETSC_ERR_SUP,"Only support for nonlinear problems");
765:   }
766:   ts->setup           = TSSetUp_PVode_Nonlinear;
767:   ts->step            = TSStep_PVode_Nonlinear;
768:   ts->setfromoptions  = TSSetFromOptions_PVode_Nonlinear;

770:   ierr  = PetscNew(TS_PVode,&cvode);
771:   ierr  = PetscMemzero(cvode,sizeof(TS_PVode));
772:   ierr  = PCCreate(ts->comm,&cvode->pc);
773:   PetscLogObjectParent(ts,cvode->pc);
774:   ts->data          = (void*)cvode;
775:   cvode->cvode_type = BDF;
776:   cvode->gtype      = PVODE_UNMODIFIED_GS;
777:   cvode->restart    = 5;
778:   cvode->linear_tol = .05;

780:   cvode->exact_final_time = PETSC_TRUE;

782:   /* set tolerance for PVode */
783:   cvode->abstol = 1e-6;
784:   cvode->reltol = 1e-6;

786:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetType_C","TSPVodeSetType_PVode",
787:                     TSPVodeSetType_PVode);
788:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGMRESRestart_C",
789:                     "TSPVodeSetGMRESRestart_PVode",
790:                     TSPVodeSetGMRESRestart_PVode);
791:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetLinearTolerance_C",
792:                     "TSPVodeSetLinearTolerance_PVode",
793:                      TSPVodeSetLinearTolerance_PVode);
794:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",
795:                     "TSPVodeSetGramSchmidtType_PVode",
796:                      TSPVodeSetGramSchmidtType_PVode);
797:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetTolerance_C",
798:                     "TSPVodeSetTolerance_PVode",
799:                      TSPVodeSetTolerance_PVode);
800:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetPC_C",
801:                     "TSPVodeGetPC_PVode",
802:                      TSPVodeGetPC_PVode);
803:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetIterations_C",
804:                     "TSPVodeGetIterations_PVode",
805:                      TSPVodeGetIterations_PVode);
806:   PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetExactFinalTime_C",
807:                     "TSPVodeSetExactFinalTime_PVode",
808:                      TSPVodeSetExactFinalTime_PVode);
809:   return(0);
810: }
811: EXTERN_C_END