Actual source code: petscpvode.c
1: /*$Id: petscpvode.c,v 1.71 2001/08/07 03:04:23 balay Exp $*/
3: /*
4: Provides a PETSc interface to PVODE. Alan Hindmarsh's parallel ODE
5: solver.
6: */
8: #include "src/ts/impls/implicit/pvode/petscpvode.h" /*I "petscts.h" I*/
10: /*
11: TSPrecond_PVode - function that we provide to PVODE to
12: evaluate the preconditioner.
14: Contributed by: Liyang Xu
16: */
17: #undef __FUNCT__
19: int TSPrecond_PVode(integertype N,realtype tn,N_Vector y,N_Vector fy,booleantype jok,
20: booleantype *jcurPtr,realtype _gamma,N_Vector ewt,realtype h,
21: realtype uround,long int *nfePtr,void *P_data,
22: N_Vector vtemp1,N_Vector vtemp2,N_Vector vtemp3)
23: {
24: TS ts = (TS) P_data;
25: TS_PVode *cvode = (TS_PVode*)ts->data;
26: PC pc = cvode->pc;
27: int ierr;
28: Mat Jac = ts->B;
29: Vec tmpy = cvode->w1;
30: PetscScalar one = 1.0,gm;
31: MatStructure str = DIFFERENT_NONZERO_PATTERN;
32:
34: /* This allows us to construct preconditioners in-place if we like */
35: MatSetUnfactored(Jac);
37: /*
38: jok - TRUE means reuse current Jacobian else recompute Jacobian
39: */
40: if (jok) {
41: ierr = MatCopy(cvode->pmat,Jac,SAME_NONZERO_PATTERN);
42: str = SAME_NONZERO_PATTERN;
43: *jcurPtr = FALSE;
44: } else {
45: /* make PETSc vector tmpy point to PVODE vector y */
46: VecPlaceArray(tmpy,N_VGetData(y));
48: /* compute the Jacobian */
49: TSComputeRHSJacobian(ts,ts->ptime,tmpy,&Jac,&Jac,&str);
51: /* copy the Jacobian matrix */
52: if (!cvode->pmat) {
53: MatDuplicate(Jac,MAT_COPY_VALUES,&cvode->pmat);
54: PetscLogObjectParent(ts,cvode->pmat);
55: }
56: MatCopy(Jac,cvode->pmat,SAME_NONZERO_PATTERN);
58: *jcurPtr = TRUE;
59: }
61: /* construct I-gamma*Jac */
62: gm = -_gamma;
63: MatScale(&gm,Jac);
64: MatShift(&one,Jac);
65:
66: PCSetOperators(pc,Jac,Jac,str);
67: return(0);
68: }
70: /*
71: TSPSolve_PVode - routine that we provide to PVode that applies the preconditioner.
72:
73: Contributed by: Liyang Xu
75: */
76: #undef __FUNCT__
78: int TSPSolve_PVode(integertype N,realtype tn,N_Vector y,N_Vector fy,N_Vector vtemp,
79: realtype _gamma,N_Vector ewt,realtype delta,long int *nfePtr,
80: N_Vector r,int lr,void *P_data,N_Vector z)
81: {
82: TS ts = (TS) P_data;
83: TS_PVode *cvode = (TS_PVode*)ts->data;
84: PC pc = cvode->pc;
85: Vec rr = cvode->w1,xx = cvode->w2;
86: int ierr;
89: /*
90: Make the PETSc work vectors rr and xx point to the arrays in the PVODE vectors
91: */
92: VecPlaceArray(rr,N_VGetData(r));
93: VecPlaceArray(xx,N_VGetData(z));
95: /*
96: Solve the Px=r and put the result in xx
97: */
98: PCApply(pc,rr,xx);
99: cvode->linear_solves++;
102: return(0);
103: }
105: /*
106: TSFunction_PVode - routine that we provide to PVode that applies the right hand side.
107:
108: Contributed by: Liyang Xu
109: */
110: #undef __FUNCT__
112: void TSFunction_PVode(int N,double t,N_Vector y,N_Vector ydot,void *ctx)
113: {
114: TS ts = (TS) ctx;
115: TS_PVode *cvode = (TS_PVode*)ts->data;
116: Vec tmpx = cvode->w1,tmpy = cvode->w2;
117: int ierr;
120: /*
121: Make the PETSc work vectors tmpx and tmpy point to the arrays in the PVODE vectors
122: */
123: VecPlaceArray(tmpx,N_VGetData(y));
124: if (ierr) {
125: (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
126: }
127: VecPlaceArray(tmpy,N_VGetData(ydot));
128: if (ierr) {
129: (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
130: }
132: /* now compute the right hand side function */
133: TSComputeRHSFunction(ts,t,tmpx,tmpy);
134: if (ierr) {
135: (*PetscErrorPrintf)("TSFunction_PVode:Could not compute RHS function. Error code %d",ierr);
136: }
137: }
139: /*
140: TSStep_PVode_Nonlinear - Calls PVode to integrate the ODE.
142: Contributed by: Liyang Xu
143: */
144: #undef __FUNCT__
146: /*
147: TSStep_PVode_Nonlinear -
148:
149: steps - number of time steps
150: time - time that integrater is terminated.
152: */
153: int TSStep_PVode_Nonlinear(TS ts,int *steps,double *time)
154: {
155: TS_PVode *cvode = (TS_PVode*)ts->data;
156: Vec sol = ts->vec_sol;
157: int ierr,i,max_steps = ts->max_steps,flag;
158: double t,tout;
159: realtype *tmp;
162: /* initialize the number of steps */
163: *steps = -ts->steps;
164: ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);
166: /* call CVSpgmr to use GMRES as the linear solver. */
167: /* setup the ode integrator with the given preconditioner */
168: CVSpgmr(cvode->mem,LEFT,cvode->gtype,cvode->restart,cvode->linear_tol,TSPrecond_PVode,TSPSolve_PVode,ts,0,0);
170: tout = ts->max_time;
171: for (i=0; i<max_steps; i++) {
172: if (ts->ptime >= tout) break;
173: VecGetArray(ts->vec_sol,&tmp);
174: N_VSetData(tmp,cvode->y);
175: flag = CVode(cvode->mem,tout,cvode->y,&t,ONE_STEP);
176: cvode->nonlinear_solves += cvode->iopt[NNI];
177: VecRestoreArray(ts->vec_sol,PETSC_NULL);
178: if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE failed");
180: if (t > tout && cvode->exact_final_time) {
181: /* interpolate to final requested time */
182: VecGetArray(ts->vec_sol,&tmp);
183: N_VSetData(tmp,cvode->y);
184: flag = CVodeDky(cvode->mem,tout,0,cvode->y);
185: VecRestoreArray(ts->vec_sol,PETSC_NULL);
186: if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE interpolation to final time failed");
187: t = tout;
188: }
190: ts->time_step = t - ts->ptime;
191: ts->ptime = t;
193: /*
194: copy the solution from cvode->y to cvode->update and sol
195: */
196: VecPlaceArray(cvode->w1,N_VGetData(cvode->y));
197: VecCopy(cvode->w1,cvode->update);
198: VecCopy(cvode->update,sol);
199:
200: ts->steps++;
201: TSMonitor(ts,ts->steps,t,sol);
202: ts->nonlinear_its = cvode->iopt[NNI];
203: ts->linear_its = cvode->iopt[SPGMR_NLI];
204: }
206: *steps += ts->steps;
207: *time = t;
209: return(0);
210: }
212: /*
214: Contributed by: Liyang Xu
215: */
216: #undef __FUNCT__
218: int TSDestroy_PVode(TS ts)
219: {
220: TS_PVode *cvode = (TS_PVode*)ts->data;
221: int ierr;
224: if (cvode->pmat) {MatDestroy(cvode->pmat);}
225: if (cvode->pc) {PCDestroy(cvode->pc);}
226: if (cvode->update) {VecDestroy(cvode->update);}
227: if (cvode->func) {VecDestroy(cvode->func);}
228: if (cvode->rhs) {VecDestroy(cvode->rhs);}
229: if (cvode->w1) {VecDestroy(cvode->w1);}
230: if (cvode->w2) {VecDestroy(cvode->w2);}
231: PetscFree(cvode);
232: return(0);
233: }
235: /*
237: Contributed by: Liyang Xu
238: */
239: #undef __FUNCT__
241: int TSSetUp_PVode_Nonlinear(TS ts)
242: {
243: TS_PVode *cvode = (TS_PVode*)ts->data;
244: int ierr,M,locsize;
245: M_Env machEnv;
246: realtype *tmp;
249: PCSetFromOptions(cvode->pc);
250: /* get the vector size */
251: VecGetSize(ts->vec_sol,&M);
252: VecGetLocalSize(ts->vec_sol,&locsize);
254: /* allocate the memory for machEnv */
255: /* machEnv = PVInitMPI(ts->comm,locsize,M); */
256: machEnv = M_EnvInit_Parallel(ts->comm, locsize, M, 0, 0);
259: /* allocate the memory for N_Vec y */
260: cvode->y = N_VNew(M,machEnv);
261: VecGetArray(ts->vec_sol,&tmp);
262: N_VSetData(tmp,cvode->y);
263: VecRestoreArray(ts->vec_sol,PETSC_NULL);
265: /* initializing vector update and func */
266: VecDuplicate(ts->vec_sol,&cvode->update);
267: VecDuplicate(ts->vec_sol,&cvode->func);
268: PetscLogObjectParent(ts,cvode->update);
269: PetscLogObjectParent(ts,cvode->func);
271: /*
272: Create work vectors for the TSPSolve_PVode() routine. Note these are
273: allocated with zero space arrays because the actual array space is provided
274: by PVode and set using VecPlaceArray().
275: */
276: VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w1);
277: VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w2);
278: PetscLogObjectParent(ts,cvode->w1);
279: PetscLogObjectParent(ts,cvode->w2);
281: PCSetVector(cvode->pc,ts->vec_sol);
283: /* allocate memory for PVode */
284: VecGetArray(ts->vec_sol,&tmp);
285: N_VSetData(tmp,cvode->y);
286: cvode->mem = CVodeMalloc(M,TSFunction_PVode,ts->ptime,cvode->y,cvode->cvode_type,
287: NEWTON,SS,&cvode->reltol,&cvode->abstol,ts,NULL,FALSE,cvode->iopt,
288: cvode->ropt,machEnv);
289: VecRestoreArray(ts->vec_sol,PETSC_NULL);
290: return(0);
291: }
293: /*
295: Contributed by: Liyang Xu
296: */
297: #undef __FUNCT__
299: int TSSetFromOptions_PVode_Nonlinear(TS ts)
300: {
301: TS_PVode *cvode = (TS_PVode*)ts->data;
302: int ierr;
303: char method[128],*btype[] = {"bdf","adams"},*otype[] = {"modified","unmodified"};
304: PetscTruth flag;
307: PetscOptionsHead("PVODE ODE solver options");
308: PetscOptionsEList("-ts_pvode_type","Scheme","TSPVodeSetType",btype,2,"bdf",method,127,&flag);
309: if (flag) {
310: PetscTruth isbdf,isadams;
312: PetscStrcmp(method,btype[0],&isbdf);
313: PetscStrcmp(method,btype[1],&isadams);
314: if (isbdf) {
315: TSPVodeSetType(ts,PVODE_BDF);
316: } else if (isadams) {
317: TSPVodeSetType(ts,PVODE_ADAMS);
318: } else {
319: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode type.n");
320: }
321: }
322: PetscOptionsEList("-ts_pvode_gramschmidt_type","Type of orthogonalization","TSPVodeSetGramSchmidtType",otype,2,"unmodified",method,127,&flag);
323: if (flag) {
324: PetscTruth ismodified,isunmodified;
326: PetscStrcmp(method,otype[0],&ismodified);
327: PetscStrcmp(method,otype[1],&isunmodified);
328: if (ismodified) {
329: TSPVodeSetGramSchmidtType(ts,PVODE_MODIFIED_GS);
330: } else if (isunmodified) {
331: TSPVodeSetGramSchmidtType(ts,PVODE_UNMODIFIED_GS);
332: } else {
333: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode Gram-Schmidt orthogonalization type n");
334: }
335: }
336: PetscOptionsReal("-ts_pvode_atol","Absolute tolerance for convergence","TSPVodeSetTolerance",cvode->abstol,&cvode->abstol,PETSC_NULL);
337: PetscOptionsReal("-ts_pvode_rtol","Relative tolerance for convergence","TSPVodeSetTolerance",cvode->reltol,&cvode->reltol,PETSC_NULL);
338: PetscOptionsReal("-ts_pvode_linear_tolerance","Convergence tolerance for linear solve","TSPVodeSetLinearTolerance",cvode->linear_tol,&cvode->linear_tol,&flag);
339: PetscOptionsInt("-ts_pvode_gmres_restart","Number of GMRES orthogonalization directions","TSPVodeSetGMRESRestart",cvode->restart,&cvode->restart,&flag);
340: PetscOptionsName("-ts_pvode_not_exact_final_time","Allow PVODE to stop near the final time, not exactly on it","TSPVodeSetExactFinalTime",&cvode->exact_final_time);
341: PetscOptionsTail();
343: return(0);
344: }
346: /*
348: Contributed by: Liyang Xu
349: */
350: #undef __FUNCT__
352: int TSPrintHelp_PVode(TS ts,char *p)
353: {
354: int ierr;
357: (*PetscHelpPrintf)(ts->comm," Options for TSPVODE integrater:n");
358: (*PetscHelpPrintf)(ts->comm," -ts_pvode_type <bdf,adams>: integration approachn",p);
359: (*PetscHelpPrintf)(ts->comm," -ts_pvode_atol aabs: absolute tolerance of ODE solutionn",p);
360: (*PetscHelpPrintf)(ts->comm," -ts_pvode_rtol rel: relative tolerance of ODE solutionn",p);
361: (*PetscHelpPrintf)(ts->comm," -ts_pvode_gramschmidt_type <unmodified,modified>n");
362: (*PetscHelpPrintf)(ts->comm," -ts_pvode_gmres_restart <restart_size> (also max. GMRES its)n");
363: (*PetscHelpPrintf)(ts->comm," -ts_pvode_linear_tolerance <tol>n");
364: (*PetscHelpPrintf)(ts->comm," -ts_pvode_not_exact_final_timen");
366: return(0);
367: }
369: /*
371: Contributed by: Liyang Xu
372: */
373: #undef __FUNCT__
375: int TSView_PVode(TS ts,PetscViewer viewer)
376: {
377: TS_PVode *cvode = (TS_PVode*)ts->data;
378: int ierr;
379: char *type;
380: PetscTruth isascii,isstring;
383: if (cvode->cvode_type == PVODE_ADAMS) {type = "Adams";}
384: else {type = "BDF: backward differentiation formula";}
386: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
387: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
388: if (isascii) {
389: PetscViewerASCIIPrintf(viewer,"PVode integrater does not use SNES!n");
390: PetscViewerASCIIPrintf(viewer,"PVode integrater type %sn",type);
391: PetscViewerASCIIPrintf(viewer,"PVode abs tol %g rel tol %gn",cvode->abstol,cvode->reltol);
392: PetscViewerASCIIPrintf(viewer,"PVode linear solver tolerance factor %gn",cvode->linear_tol);
393: PetscViewerASCIIPrintf(viewer,"PVode GMRES max iterations (same as restart in PVODE) %dn",cvode->restart);
394: if (cvode->gtype == PVODE_MODIFIED_GS) {
395: PetscViewerASCIIPrintf(viewer,"PVode using modified Gram-Schmidt for orthogonalization in GMRESn");
396: } else {
397: PetscViewerASCIIPrintf(viewer,"PVode using unmodified (classical) Gram-Schmidt for orthogonalization in GMRESn");
398: }
399: } else if (isstring) {
400: PetscViewerStringSPrintf(viewer,"Pvode type %s",type);
401: } else {
402: SETERRQ1(1,"Viewer type %s not supported by TS PVode",((PetscObject)viewer)->type_name);
403: }
404: PetscViewerASCIIPushTab(viewer);
405: PCView(cvode->pc,viewer);
406: PetscViewerASCIIPopTab(viewer);
408: return(0);
409: }
412: /* --------------------------------------------------------------------------*/
413: EXTERN_C_BEGIN
414: #undef __FUNCT__
416: int TSPVodeSetType_PVode(TS ts,TSPVodeType type)
417: {
418: TS_PVode *cvode = (TS_PVode*)ts->data;
419:
421: cvode->cvode_type = type;
422: return(0);
423: }
424: EXTERN_C_END
426: EXTERN_C_BEGIN
427: #undef __FUNCT__
429: int TSPVodeSetGMRESRestart_PVode(TS ts,int restart)
430: {
431: TS_PVode *cvode = (TS_PVode*)ts->data;
432:
434: cvode->restart = restart;
435: return(0);
436: }
437: EXTERN_C_END
439: EXTERN_C_BEGIN
440: #undef __FUNCT__
442: int TSPVodeSetLinearTolerance_PVode(TS ts,double tol)
443: {
444: TS_PVode *cvode = (TS_PVode*)ts->data;
445:
447: cvode->linear_tol = tol;
448: return(0);
449: }
450: EXTERN_C_END
452: EXTERN_C_BEGIN
453: #undef __FUNCT__
455: int TSPVodeSetGramSchmidtType_PVode(TS ts,TSPVodeGramSchmidtType type)
456: {
457: TS_PVode *cvode = (TS_PVode*)ts->data;
458:
460: cvode->gtype = type;
462: return(0);
463: }
464: EXTERN_C_END
466: EXTERN_C_BEGIN
467: #undef __FUNCT__
469: int TSPVodeSetTolerance_PVode(TS ts,double aabs,double rel)
470: {
471: TS_PVode *cvode = (TS_PVode*)ts->data;
472:
474: if (aabs != PETSC_DECIDE) cvode->abstol = aabs;
475: if (rel != PETSC_DECIDE) cvode->reltol = rel;
476: return(0);
477: }
478: EXTERN_C_END
480: EXTERN_C_BEGIN
481: #undef __FUNCT__
483: int TSPVodeGetPC_PVode(TS ts,PC *pc)
484: {
485: TS_PVode *cvode = (TS_PVode*)ts->data;
488: *pc = cvode->pc;
490: return(0);
491: }
492: EXTERN_C_END
494: EXTERN_C_BEGIN
495: #undef __FUNCT__
497: int TSPVodeGetIterations_PVode(TS ts,int *nonlin,int *lin)
498: {
499: TS_PVode *cvode = (TS_PVode*)ts->data;
500:
502: if (nonlin) *nonlin = cvode->nonlinear_solves;
503: if (lin) *lin = cvode->linear_solves;
504: return(0);
505: }
506: EXTERN_C_END
507:
508: EXTERN_C_BEGIN
509: #undef __FUNCT__
511: int TSPVodeSetExactFinalTime_PVode(TS ts,PetscTruth s)
512: {
513: TS_PVode *cvode = (TS_PVode*)ts->data;
514:
516: cvode->exact_final_time = s;
517: return(0);
518: }
519: EXTERN_C_END
520: /* -------------------------------------------------------------------------------------------*/
522: #undef __FUNCT__
524: /*@C
525: TSPVodeGetIterations - Gets the number of nonlinear and linear iterations used so far by PVode.
527: Not Collective
529: Input parameters:
530: . ts - the time-step context
532: Output Parameters:
533: + nonlin - number of nonlinear iterations
534: - lin - number of linear iterations
536: Level: advanced
538: Notes:
539: These return the number since the creation of the TS object
541: .keywords: non-linear iterations, linear iterations
543: .seealso: TSPVodeSetType(), TSPVodeSetGMRESRestart(),
544: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
545: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
546: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
547: TSPVodeSetExactFinalTime()
549: @*/
550: int TSPVodeGetIterations(TS ts,int *nonlin,int *lin)
551: {
552: int ierr,(*f)(TS,int*,int*);
553:
555: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetIterations_C",(void (**)(void))&f);
556: if (f) {
557: (*f)(ts,nonlin,lin);
558: }
559: return(0);
560: }
562: #undef __FUNCT__
564: /*@
565: TSPVodeSetType - Sets the method that PVode will use for integration.
567: Collective on TS
569: Input parameters:
570: + ts - the time-step context
571: - type - one of PVODE_ADAMS or PVODE_BDF
573: Contributed by: Liyang Xu
575: Level: intermediate
577: .keywords: Adams, backward differentiation formula
579: .seealso: TSPVodeGetIterations(), TSPVodeSetGMRESRestart(),
580: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
581: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
582: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
583: TSPVodeSetExactFinalTime()
584: @*/
585: int TSPVodeSetType(TS ts,TSPVodeType type)
586: {
587: int ierr,(*f)(TS,TSPVodeType);
588:
590: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetType_C",(void (**)(void))&f);
591: if (f) {
592: (*f)(ts,type);
593: }
594: return(0);
595: }
597: #undef __FUNCT__
599: /*@
600: TSPVodeSetGMRESRestart - Sets the dimension of the Krylov space used by
601: GMRES in the linear solver in PVODE. PVODE DOES NOT use restarted GMRES so
602: this is ALSO the maximum number of GMRES steps that will be used.
604: Collective on TS
606: Input parameters:
607: + ts - the time-step context
608: - restart - number of direction vectors (the restart size).
610: Level: advanced
612: .keywords: GMRES, restart
614: .seealso: TSPVodeGetIterations(), TSPVodeSetType(),
615: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
616: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
617: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
618: TSPVodeSetExactFinalTime()
620: @*/
621: int TSPVodeSetGMRESRestart(TS ts,int restart)
622: {
623: int ierr,(*f)(TS,int);
626: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGMRESRestart_C",(void (**)(void))&f);
627: if (f) {
628: (*f)(ts,restart);
629: }
631: return(0);
632: }
634: #undef __FUNCT__
636: /*@
637: TSPVodeSetLinearTolerance - Sets the tolerance used to solve the linear
638: system by PVODE.
640: Collective on TS
642: Input parameters:
643: + ts - the time-step context
644: - tol - the factor by which the tolerance on the nonlinear solver is
645: multiplied to get the tolerance on the linear solver, .05 by default.
647: Level: advanced
649: .keywords: GMRES, linear convergence tolerance, PVODE
651: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
652: TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
653: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
654: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
655: TSPVodeSetExactFinalTime()
657: @*/
658: int TSPVodeSetLinearTolerance(TS ts,double tol)
659: {
660: int ierr,(*f)(TS,double);
661:
663: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetLinearTolerance_C",(void (**)(void))&f);
664: if (f) {
665: (*f)(ts,tol);
666: }
667: return(0);
668: }
670: #undef __FUNCT__
672: /*@
673: TSPVodeSetGramSchmidtType - Sets type of orthogonalization used
674: in GMRES method by PVODE linear solver.
676: Collective on TS
678: Input parameters:
679: + ts - the time-step context
680: - type - either PVODE_MODIFIED_GS or PVODE_CLASSICAL_GS
682: Level: advanced
684: .keywords: PVode, orthogonalization
686: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
687: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(),
688: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
689: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
690: TSPVodeSetExactFinalTime()
692: @*/
693: int TSPVodeSetGramSchmidtType(TS ts,TSPVodeGramSchmidtType type)
694: {
695: int ierr,(*f)(TS,TSPVodeGramSchmidtType);
696:
698: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",(void (**)(void))&f);
699: if (f) {
700: (*f)(ts,type);
701: }
702: return(0);
703: }
705: #undef __FUNCT__
707: /*@
708: TSPVodeSetTolerance - Sets the absolute and relative tolerance used by
709: PVode for error control.
711: Collective on TS
713: Input parameters:
714: + ts - the time-step context
715: . aabs - the absolute tolerance
716: - rel - the relative tolerance
718: Contributed by: Liyang Xu
720: Level: intermediate
722: .keywords: PVode, tolerance
724: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
725: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(),
726: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
727: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
728: TSPVodeSetExactFinalTime()
730: @*/
731: int TSPVodeSetTolerance(TS ts,double aabs,double rel)
732: {
733: int ierr,(*f)(TS,double,double);
734:
736: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetTolerance_C",(void (**)(void))&f);
737: if (f) {
738: (*f)(ts,aabs,rel);
739: }
740: return(0);
741: }
743: #undef __FUNCT__
745: /*@
746: TSPVodeGetPC - Extract the PC context from a time-step context for PVode.
748: Input Parameter:
749: . ts - the time-step context
751: Output Parameter:
752: . pc - the preconditioner context
754: Level: advanced
756: Contributed by: Liyang Xu
758: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
759: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
760: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
761: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance()
762: @*/
763: int TSPVodeGetPC(TS ts,PC *pc)
764: {
765: int ierr,(*f)(TS,PC *);
768: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetPC_C",(void (**)(void))&f);
769: if (f) {
770: (*f)(ts,pc);
771: } else {
772: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"TS must be of PVode type to extract the PC");
773: }
775: return(0);
776: }
778: #undef __FUNCT__
780: /*@
781: TSPVodeSetExactFinalTime - Determines if PVode interpolates solution to the
782: exact final time requested by the user or just returns it at the final time
783: it computed. (Defaults to true).
785: Input Parameter:
786: + ts - the time-step context
787: - ft - PETSC_TRUE if interpolates, else PETSC_FALSE
789: Level: beginner
791: .seealso:TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
792: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
793: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
794: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC()
795: @*/
796: int TSPVodeSetExactFinalTime(TS ts,PetscTruth ft)
797: {
798: int ierr,(*f)(TS,PetscTruth);
801: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetExactFinalTime_C",(void (**)(void))&f);
802: if (f) {
803: (*f)(ts,ft);
804: }
806: return(0);
807: }
809: /* -------------------------------------------------------------------------------------------*/
811: /*
813: Contributed by: Liyang Xu
814: */
815: EXTERN_C_BEGIN
816: #undef __FUNCT__
818: int TSCreate_PVode(TS ts)
819: {
820: TS_PVode *cvode;
821: int ierr;
824: ts->ops->destroy = TSDestroy_PVode;
825: ts->ops->view = TSView_PVode;
827: if (ts->problem_type != TS_NONLINEAR) {
828: SETERRQ(PETSC_ERR_SUP,"Only support for nonlinear problems");
829: }
830: ts->ops->setup = TSSetUp_PVode_Nonlinear;
831: ts->ops->step = TSStep_PVode_Nonlinear;
832: ts->ops->setfromoptions = TSSetFromOptions_PVode_Nonlinear;
834: ierr = PetscNew(TS_PVode,&cvode);
835: ierr = PetscMemzero(cvode,sizeof(TS_PVode));
836: ierr = PCCreate(ts->comm,&cvode->pc);
837: PetscLogObjectParent(ts,cvode->pc);
838: ts->data = (void*)cvode;
839: cvode->cvode_type = BDF;
840: cvode->gtype = PVODE_UNMODIFIED_GS;
841: cvode->restart = 5;
842: cvode->linear_tol = .05;
844: cvode->exact_final_time = PETSC_TRUE;
846: /* set tolerance for PVode */
847: cvode->abstol = 1e-6;
848: cvode->reltol = 1e-6;
850: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetType_C","TSPVodeSetType_PVode",
851: TSPVodeSetType_PVode);
852: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGMRESRestart_C",
853: "TSPVodeSetGMRESRestart_PVode",
854: TSPVodeSetGMRESRestart_PVode);
855: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetLinearTolerance_C",
856: "TSPVodeSetLinearTolerance_PVode",
857: TSPVodeSetLinearTolerance_PVode);
858: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",
859: "TSPVodeSetGramSchmidtType_PVode",
860: TSPVodeSetGramSchmidtType_PVode);
861: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetTolerance_C",
862: "TSPVodeSetTolerance_PVode",
863: TSPVodeSetTolerance_PVode);
864: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetPC_C",
865: "TSPVodeGetPC_PVode",
866: TSPVodeGetPC_PVode);
867: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetIterations_C",
868: "TSPVodeGetIterations_PVode",
869: TSPVodeGetIterations_PVode);
870: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetExactFinalTime_C",
871: "TSPVodeSetExactFinalTime_PVode",
872: TSPVodeSetExactFinalTime_PVode);
873: return(0);
874: }
875: EXTERN_C_END