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*/
12: /*@C
13: TSPVodeGetParameters - Extract "iopt" and "ropt" PVODE parameters
15: Input Parameter:
16: . ts - the time-step context
18: Output Parameters:
19: + opt_size - size of the parameter arrays (will be set to PVODE value OPT_SIZE)
20: . iopt - the integer parameters
21: - ropt - the double paramters
24: Level: advanced
25: g
26: Notes: You may pass PETSC_NULL for any value you do not desire
28: PETSc initializes these array with the default PVODE values of 0, you may
29: change them before calling the TS solver.
31: See the PVODE include file cvode.h for the definitions of the fields
33: Suggested by: Timothy William Chevalier
35: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
36: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
37: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
38: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance()
39: @*/
40: int TSPVodeGetParameters(TS ts,int *opt_size,long int *iopt[],double *ropt[])
41: {
42: TS_PVode *cvode = (TS_PVode*)ts->data;
45: if (opt_size) *opt_size = OPT_SIZE;
46: if (iopt) *iopt = cvode->iopt;
47: if (ropt) *ropt = cvode->ropt;
48: return(0);
49: }
51: /*
52: TSPrecond_PVode - function that we provide to PVODE to
53: evaluate the preconditioner.
55: Contributed by: Liyang Xu
57: */
60: int TSPrecond_PVode(integertype N,realtype tn,N_Vector y,N_Vector fy,booleantype jok,
61: booleantype *jcurPtr,realtype _gamma,N_Vector ewt,realtype h,
62: realtype uround,long int *nfePtr,void *P_data,
63: N_Vector vtemp1,N_Vector vtemp2,N_Vector vtemp3)
64: {
65: TS ts = (TS) P_data;
66: TS_PVode *cvode = (TS_PVode*)ts->data;
67: PC pc = cvode->pc;
68: int ierr;
69: Mat Jac = ts->B;
70: Vec tmpy = cvode->w1;
71: PetscScalar one = 1.0,gm;
72: MatStructure str = DIFFERENT_NONZERO_PATTERN;
73:
75: /* This allows us to construct preconditioners in-place if we like */
76: MatSetUnfactored(Jac);
78: /*
79: jok - TRUE means reuse current Jacobian else recompute Jacobian
80: */
81: if (jok) {
82: MatCopy(cvode->pmat,Jac,SAME_NONZERO_PATTERN);
83: str = SAME_NONZERO_PATTERN;
84: *jcurPtr = FALSE;
85: } else {
86: /* make PETSc vector tmpy point to PVODE vector y */
87: VecPlaceArray(tmpy,N_VGetData(y));
89: /* compute the Jacobian */
90: TSComputeRHSJacobian(ts,ts->ptime,tmpy,&Jac,&Jac,&str);
92: /* copy the Jacobian matrix */
93: if (!cvode->pmat) {
94: MatDuplicate(Jac,MAT_COPY_VALUES,&cvode->pmat);
95: PetscLogObjectParent(ts,cvode->pmat);
96: }
97: MatCopy(Jac,cvode->pmat,SAME_NONZERO_PATTERN);
99: *jcurPtr = TRUE;
100: }
102: /* construct I-gamma*Jac */
103: gm = -_gamma;
104: MatScale(&gm,Jac);
105: MatShift(&one,Jac);
106:
107: PCSetOperators(pc,Jac,Jac,str);
108: return(0);
109: }
111: /*
112: TSPSolve_PVode - routine that we provide to PVode that applies the preconditioner.
113:
114: Contributed by: Liyang Xu
116: */
119: int TSPSolve_PVode(integertype N,realtype tn,N_Vector y,N_Vector fy,N_Vector vtemp,
120: realtype _gamma,N_Vector ewt,realtype delta,long int *nfePtr,
121: N_Vector r,int lr,void *P_data,N_Vector z)
122: {
123: TS ts = (TS) P_data;
124: TS_PVode *cvode = (TS_PVode*)ts->data;
125: PC pc = cvode->pc;
126: Vec rr = cvode->w1,xx = cvode->w2;
127: int ierr;
130: /*
131: Make the PETSc work vectors rr and xx point to the arrays in the PVODE vectors
132: */
133: VecPlaceArray(rr,N_VGetData(r));
134: VecPlaceArray(xx,N_VGetData(z));
136: /*
137: Solve the Px=r and put the result in xx
138: */
139: PCApply(pc,rr,xx,PC_LEFT);
140: cvode->linear_solves++;
143: return(0);
144: }
146: /*
147: TSFunction_PVode - routine that we provide to PVode that applies the right hand side.
148:
149: Contributed by: Liyang Xu
150: */
153: void TSFunction_PVode(int N,double t,N_Vector y,N_Vector ydot,void *ctx)
154: {
155: TS ts = (TS) ctx;
156: TS_PVode *cvode = (TS_PVode*)ts->data;
157: Vec tmpx = cvode->w1,tmpy = cvode->w2;
158: int ierr;
161: /*
162: Make the PETSc work vectors tmpx and tmpy point to the arrays in the PVODE vectors
163: */
164: VecPlaceArray(tmpx,N_VGetData(y));
165: if (ierr) {
166: (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
167: }
168: VecPlaceArray(tmpy,N_VGetData(ydot));
169: if (ierr) {
170: (*PetscErrorPrintf)("TSFunction_PVode:Could not place array. Error code %d",ierr);
171: }
173: /* now compute the right hand side function */
174: TSComputeRHSFunction(ts,t,tmpx,tmpy);
175: if (ierr) {
176: (*PetscErrorPrintf)("TSFunction_PVode:Could not compute RHS function. Error code %d",ierr);
177: }
178: }
180: /*
181: TSStep_PVode_Nonlinear - Calls PVode to integrate the ODE.
183: Contributed by: Liyang Xu
184: */
187: /*
188: TSStep_PVode_Nonlinear -
189:
190: steps - number of time steps
191: time - time that integrater is terminated.
193: */
194: int TSStep_PVode_Nonlinear(TS ts,int *steps,double *time)
195: {
196: TS_PVode *cvode = (TS_PVode*)ts->data;
197: Vec sol = ts->vec_sol;
198: int ierr,i,max_steps = ts->max_steps,flag;
199: double t,tout;
200: realtype *tmp;
203: /* initialize the number of steps */
204: *steps = -ts->steps;
205: TSMonitor(ts,ts->steps,ts->ptime,sol);
207: /* call CVSpgmr to use GMRES as the linear solver. */
208: /* setup the ode integrator with the given preconditioner */
209: CVSpgmr(cvode->mem,LEFT,cvode->gtype,cvode->restart,cvode->linear_tol,TSPrecond_PVode,TSPSolve_PVode,ts,0,0);
211: tout = ts->max_time;
212: for (i=0; i<max_steps; i++) {
213: if (ts->ptime >= tout) break;
214: VecGetArray(ts->vec_sol,&tmp);
215: N_VSetData(tmp,cvode->y);
216: flag = CVode(cvode->mem,tout,cvode->y,&t,ONE_STEP);
217: cvode->nonlinear_solves += cvode->iopt[NNI];
218: VecRestoreArray(ts->vec_sol,PETSC_NULL);
219: if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE failed");
221: if (t > tout && cvode->exact_final_time) {
222: /* interpolate to final requested time */
223: VecGetArray(ts->vec_sol,&tmp);
224: N_VSetData(tmp,cvode->y);
225: flag = CVodeDky(cvode->mem,tout,0,cvode->y);
226: VecRestoreArray(ts->vec_sol,PETSC_NULL);
227: if (flag != SUCCESS) SETERRQ(PETSC_ERR_LIB,"PVODE interpolation to final time failed");
228: t = tout;
229: }
231: ts->time_step = t - ts->ptime;
232: ts->ptime = t;
234: /*
235: copy the solution from cvode->y to cvode->update and sol
236: */
237: VecPlaceArray(cvode->w1,N_VGetData(cvode->y));
238: VecCopy(cvode->w1,cvode->update);
239: VecCopy(cvode->update,sol);
240:
241: ts->steps++;
242: TSMonitor(ts,ts->steps,t,sol);
243: ts->nonlinear_its = cvode->iopt[NNI];
244: ts->linear_its = cvode->iopt[SPGMR_NLI];
245: }
247: *steps += ts->steps;
248: *time = t;
250: return(0);
251: }
253: /*
255: Contributed by: Liyang Xu
256: */
259: int TSDestroy_PVode(TS ts)
260: {
261: TS_PVode *cvode = (TS_PVode*)ts->data;
262: int ierr;
265: if (cvode->pmat) {MatDestroy(cvode->pmat);}
266: if (cvode->pc) {PCDestroy(cvode->pc);}
267: if (cvode->update) {VecDestroy(cvode->update);}
268: if (cvode->func) {VecDestroy(cvode->func);}
269: if (cvode->rhs) {VecDestroy(cvode->rhs);}
270: if (cvode->w1) {VecDestroy(cvode->w1);}
271: if (cvode->w2) {VecDestroy(cvode->w2);}
272: MPI_Comm_free(&(cvode->comm_pvode));
273: PetscFree(cvode);
274: return(0);
275: }
277: /*
279: Contributed by: Liyang Xu
280: */
283: int TSSetUp_PVode_Nonlinear(TS ts)
284: {
285: TS_PVode *cvode = (TS_PVode*)ts->data;
286: int ierr,M,locsize;
287: M_Env machEnv;
288: realtype *tmp;
291: PCSetFromOptions(cvode->pc);
292: /* get the vector size */
293: VecGetSize(ts->vec_sol,&M);
294: VecGetLocalSize(ts->vec_sol,&locsize);
296: /* allocate the memory for machEnv */
297: /* machEnv = PVInitMPI(cvode->>comm_pvode,locsize,M); */
298: machEnv = M_EnvInit_Parallel(cvode->comm_pvode, locsize, M, 0, 0);
301: /* allocate the memory for N_Vec y */
302: cvode->y = N_VNew(M,machEnv);
303: VecGetArray(ts->vec_sol,&tmp);
304: N_VSetData(tmp,cvode->y);
305: VecRestoreArray(ts->vec_sol,PETSC_NULL);
307: /* initializing vector update and func */
308: VecDuplicate(ts->vec_sol,&cvode->update);
309: VecDuplicate(ts->vec_sol,&cvode->func);
310: PetscLogObjectParent(ts,cvode->update);
311: PetscLogObjectParent(ts,cvode->func);
313: /*
314: Create work vectors for the TSPSolve_PVode() routine. Note these are
315: allocated with zero space arrays because the actual array space is provided
316: by PVode and set using VecPlaceArray().
317: */
318: VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w1);
319: VecCreateMPIWithArray(ts->comm,locsize,PETSC_DECIDE,0,&cvode->w2);
320: PetscLogObjectParent(ts,cvode->w1);
321: PetscLogObjectParent(ts,cvode->w2);
323: PCSetVector(cvode->pc,ts->vec_sol);
325: /* allocate memory for PVode */
326: VecGetArray(ts->vec_sol,&tmp);
327: N_VSetData(tmp,cvode->y);
328: cvode->mem = CVodeMalloc(M,TSFunction_PVode,ts->ptime,cvode->y,cvode->cvode_type,
329: NEWTON,SS,&cvode->reltol,&cvode->abstol,ts,NULL,TRUE,cvode->iopt,
330: cvode->ropt,machEnv);
331: VecRestoreArray(ts->vec_sol,PETSC_NULL);
332: return(0);
333: }
335: /*
337: Contributed by: Liyang Xu
338: */
341: int TSSetFromOptions_PVode_Nonlinear(TS ts)
342: {
343: TS_PVode *cvode = (TS_PVode*)ts->data;
344: int ierr;
345: char method[128],*btype[] = {"bdf","adams"},*otype[] = {"modified","unmodified"};
346: PetscTruth flag;
349: PetscOptionsHead("PVODE ODE solver options");
350: PetscOptionsEList("-ts_pvode_type","Scheme","TSPVodeSetType",btype,2,"bdf",method,127,&flag);
351: if (flag) {
352: PetscTruth isbdf,isadams;
354: PetscStrcmp(method,btype[0],&isbdf);
355: PetscStrcmp(method,btype[1],&isadams);
356: if (isbdf) {
357: TSPVodeSetType(ts,PVODE_BDF);
358: } else if (isadams) {
359: TSPVodeSetType(ts,PVODE_ADAMS);
360: } else {
361: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode type.\n");
362: }
363: }
364: PetscOptionsEList("-ts_pvode_gramschmidt_type","Type of orthogonalization","TSPVodeSetGramSchmidtType",otype,2,"unmodified",method,127,&flag);
365: if (flag) {
366: PetscTruth ismodified,isunmodified;
368: PetscStrcmp(method,otype[0],&ismodified);
369: PetscStrcmp(method,otype[1],&isunmodified);
370: if (ismodified) {
371: TSPVodeSetGramSchmidtType(ts,PVODE_MODIFIED_GS);
372: } else if (isunmodified) {
373: TSPVodeSetGramSchmidtType(ts,PVODE_UNMODIFIED_GS);
374: } else {
375: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown PVode Gram-Schmidt orthogonalization type \n");
376: }
377: }
378: PetscOptionsReal("-ts_pvode_atol","Absolute tolerance for convergence","TSPVodeSetTolerance",cvode->abstol,&cvode->abstol,PETSC_NULL);
379: PetscOptionsReal("-ts_pvode_rtol","Relative tolerance for convergence","TSPVodeSetTolerance",cvode->reltol,&cvode->reltol,PETSC_NULL);
380: PetscOptionsReal("-ts_pvode_linear_tolerance","Convergence tolerance for linear solve","TSPVodeSetLinearTolerance",cvode->linear_tol,&cvode->linear_tol,&flag);
381: PetscOptionsInt("-ts_pvode_gmres_restart","Number of GMRES orthogonalization directions","TSPVodeSetGMRESRestart",cvode->restart,&cvode->restart,&flag);
382: PetscOptionsName("-ts_pvode_not_exact_final_time","Allow PVODE to stop near the final time, not exactly on it","TSPVodeSetExactFinalTime",&cvode->exact_final_time);
383: PetscOptionsTail();
385: return(0);
386: }
388: /*
390: Contributed by: Liyang Xu
391: */
394: int TSPrintHelp_PVode(TS ts,char *p)
395: {
396: int ierr;
399: (*PetscHelpPrintf)(ts->comm," Options for TSPVODE integrater:\n");
400: (*PetscHelpPrintf)(ts->comm," -ts_pvode_type <bdf,adams>: integration approach\n",p);
401: (*PetscHelpPrintf)(ts->comm," -ts_pvode_atol aabs: absolute tolerance of ODE solution\n",p);
402: (*PetscHelpPrintf)(ts->comm," -ts_pvode_rtol rel: relative tolerance of ODE solution\n",p);
403: (*PetscHelpPrintf)(ts->comm," -ts_pvode_gramschmidt_type <unmodified,modified>\n");
404: (*PetscHelpPrintf)(ts->comm," -ts_pvode_gmres_restart <restart_size> (also max. GMRES its)\n");
405: (*PetscHelpPrintf)(ts->comm," -ts_pvode_linear_tolerance <tol>\n");
406: (*PetscHelpPrintf)(ts->comm," -ts_pvode_not_exact_final_time\n");
408: return(0);
409: }
411: /*
413: Contributed by: Liyang Xu
414: */
417: int TSView_PVode(TS ts,PetscViewer viewer)
418: {
419: TS_PVode *cvode = (TS_PVode*)ts->data;
420: int ierr;
421: char *type;
422: PetscTruth isascii,isstring;
425: if (cvode->cvode_type == PVODE_ADAMS) {type = "Adams";}
426: else {type = "BDF: backward differentiation formula";}
428: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
429: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
430: if (isascii) {
431: PetscViewerASCIIPrintf(viewer,"PVode integrater does not use SNES!\n");
432: PetscViewerASCIIPrintf(viewer,"PVode integrater type %s\n",type);
433: PetscViewerASCIIPrintf(viewer,"PVode abs tol %g rel tol %g\n",cvode->abstol,cvode->reltol);
434: PetscViewerASCIIPrintf(viewer,"PVode linear solver tolerance factor %g\n",cvode->linear_tol);
435: PetscViewerASCIIPrintf(viewer,"PVode GMRES max iterations (same as restart in PVODE) %d\n",cvode->restart);
436: if (cvode->gtype == PVODE_MODIFIED_GS) {
437: PetscViewerASCIIPrintf(viewer,"PVode using modified Gram-Schmidt for orthogonalization in GMRES\n");
438: } else {
439: PetscViewerASCIIPrintf(viewer,"PVode using unmodified (classical) Gram-Schmidt for orthogonalization in GMRES\n");
440: }
441: } else if (isstring) {
442: PetscViewerStringSPrintf(viewer,"Pvode type %s",type);
443: } else {
444: SETERRQ1(1,"Viewer type %s not supported by TS PVode",((PetscObject)viewer)->type_name);
445: }
446: PetscViewerASCIIPushTab(viewer);
447: PCView(cvode->pc,viewer);
448: PetscViewerASCIIPopTab(viewer);
450: return(0);
451: }
454: /* --------------------------------------------------------------------------*/
455: EXTERN_C_BEGIN
458: int TSPVodeSetType_PVode(TS ts,TSPVodeType type)
459: {
460: TS_PVode *cvode = (TS_PVode*)ts->data;
461:
463: cvode->cvode_type = type;
464: return(0);
465: }
466: EXTERN_C_END
468: EXTERN_C_BEGIN
471: int TSPVodeSetGMRESRestart_PVode(TS ts,int restart)
472: {
473: TS_PVode *cvode = (TS_PVode*)ts->data;
474:
476: cvode->restart = restart;
477: return(0);
478: }
479: EXTERN_C_END
481: EXTERN_C_BEGIN
484: int TSPVodeSetLinearTolerance_PVode(TS ts,double tol)
485: {
486: TS_PVode *cvode = (TS_PVode*)ts->data;
487:
489: cvode->linear_tol = tol;
490: return(0);
491: }
492: EXTERN_C_END
494: EXTERN_C_BEGIN
497: int TSPVodeSetGramSchmidtType_PVode(TS ts,TSPVodeGramSchmidtType type)
498: {
499: TS_PVode *cvode = (TS_PVode*)ts->data;
500:
502: cvode->gtype = type;
504: return(0);
505: }
506: EXTERN_C_END
508: EXTERN_C_BEGIN
511: int TSPVodeSetTolerance_PVode(TS ts,double aabs,double rel)
512: {
513: TS_PVode *cvode = (TS_PVode*)ts->data;
514:
516: if (aabs != PETSC_DECIDE) cvode->abstol = aabs;
517: if (rel != PETSC_DECIDE) cvode->reltol = rel;
518: return(0);
519: }
520: EXTERN_C_END
522: EXTERN_C_BEGIN
525: int TSPVodeGetPC_PVode(TS ts,PC *pc)
526: {
527: TS_PVode *cvode = (TS_PVode*)ts->data;
530: *pc = cvode->pc;
532: return(0);
533: }
534: EXTERN_C_END
536: EXTERN_C_BEGIN
539: int TSPVodeGetIterations_PVode(TS ts,int *nonlin,int *lin)
540: {
541: TS_PVode *cvode = (TS_PVode*)ts->data;
542:
544: if (nonlin) *nonlin = cvode->nonlinear_solves;
545: if (lin) *lin = cvode->linear_solves;
546: return(0);
547: }
548: EXTERN_C_END
549:
550: EXTERN_C_BEGIN
553: int TSPVodeSetExactFinalTime_PVode(TS ts,PetscTruth s)
554: {
555: TS_PVode *cvode = (TS_PVode*)ts->data;
556:
558: cvode->exact_final_time = s;
559: return(0);
560: }
561: EXTERN_C_END
562: /* -------------------------------------------------------------------------------------------*/
566: /*@C
567: TSPVodeGetIterations - Gets the number of nonlinear and linear iterations used so far by PVode.
569: Not Collective
571: Input parameters:
572: . ts - the time-step context
574: Output Parameters:
575: + nonlin - number of nonlinear iterations
576: - lin - number of linear iterations
578: Level: advanced
580: Notes:
581: These return the number since the creation of the TS object
583: .keywords: non-linear iterations, linear iterations
585: .seealso: TSPVodeSetType(), TSPVodeSetGMRESRestart(),
586: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
587: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
588: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
589: TSPVodeSetExactFinalTime()
591: @*/
592: int TSPVodeGetIterations(TS ts,int *nonlin,int *lin)
593: {
594: int ierr,(*f)(TS,int*,int*);
595:
597: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetIterations_C",(void (**)(void))&f);
598: if (f) {
599: (*f)(ts,nonlin,lin);
600: }
601: return(0);
602: }
606: /*@
607: TSPVodeSetType - Sets the method that PVode will use for integration.
609: Collective on TS
611: Input parameters:
612: + ts - the time-step context
613: - type - one of PVODE_ADAMS or PVODE_BDF
615: Contributed by: Liyang Xu
617: Level: intermediate
619: .keywords: Adams, backward differentiation formula
621: .seealso: TSPVodeGetIterations(), TSPVodeSetGMRESRestart(),
622: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
623: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
624: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
625: TSPVodeSetExactFinalTime()
626: @*/
627: int TSPVodeSetType(TS ts,TSPVodeType type)
628: {
629: int ierr,(*f)(TS,TSPVodeType);
630:
632: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetType_C",(void (**)(void))&f);
633: if (f) {
634: (*f)(ts,type);
635: }
636: return(0);
637: }
641: /*@
642: TSPVodeSetGMRESRestart - Sets the dimension of the Krylov space used by
643: GMRES in the linear solver in PVODE. PVODE DOES NOT use restarted GMRES so
644: this is ALSO the maximum number of GMRES steps that will be used.
646: Collective on TS
648: Input parameters:
649: + ts - the time-step context
650: - restart - number of direction vectors (the restart size).
652: Level: advanced
654: .keywords: GMRES, restart
656: .seealso: TSPVodeGetIterations(), TSPVodeSetType(),
657: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
658: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
659: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
660: TSPVodeSetExactFinalTime()
662: @*/
663: int TSPVodeSetGMRESRestart(TS ts,int restart)
664: {
665: int ierr,(*f)(TS,int);
668: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGMRESRestart_C",(void (**)(void))&f);
669: if (f) {
670: (*f)(ts,restart);
671: }
673: return(0);
674: }
678: /*@
679: TSPVodeSetLinearTolerance - Sets the tolerance used to solve the linear
680: system by PVODE.
682: Collective on TS
684: Input parameters:
685: + ts - the time-step context
686: - tol - the factor by which the tolerance on the nonlinear solver is
687: multiplied to get the tolerance on the linear solver, .05 by default.
689: Level: advanced
691: .keywords: GMRES, linear convergence tolerance, PVODE
693: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
694: TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
695: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
696: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
697: TSPVodeSetExactFinalTime()
699: @*/
700: int TSPVodeSetLinearTolerance(TS ts,double tol)
701: {
702: int ierr,(*f)(TS,double);
703:
705: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetLinearTolerance_C",(void (**)(void))&f);
706: if (f) {
707: (*f)(ts,tol);
708: }
709: return(0);
710: }
714: /*@
715: TSPVodeSetGramSchmidtType - Sets type of orthogonalization used
716: in GMRES method by PVODE linear solver.
718: Collective on TS
720: Input parameters:
721: + ts - the time-step context
722: - type - either PVODE_MODIFIED_GS or PVODE_CLASSICAL_GS
724: Level: advanced
726: .keywords: PVode, orthogonalization
728: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
729: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(),
730: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
731: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
732: TSPVodeSetExactFinalTime()
734: @*/
735: int TSPVodeSetGramSchmidtType(TS ts,TSPVodeGramSchmidtType type)
736: {
737: int ierr,(*f)(TS,TSPVodeGramSchmidtType);
738:
740: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",(void (**)(void))&f);
741: if (f) {
742: (*f)(ts,type);
743: }
744: return(0);
745: }
749: /*@
750: TSPVodeSetTolerance - Sets the absolute and relative tolerance used by
751: PVode for error control.
753: Collective on TS
755: Input parameters:
756: + ts - the time-step context
757: . aabs - the absolute tolerance
758: - rel - the relative tolerance
760: Contributed by: Liyang Xu
762: Level: intermediate
764: .keywords: PVode, tolerance
766: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
767: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(),
768: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
769: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC(),
770: TSPVodeSetExactFinalTime()
772: @*/
773: int TSPVodeSetTolerance(TS ts,double aabs,double rel)
774: {
775: int ierr,(*f)(TS,double,double);
776:
778: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetTolerance_C",(void (**)(void))&f);
779: if (f) {
780: (*f)(ts,aabs,rel);
781: }
782: return(0);
783: }
787: /*@
788: TSPVodeGetPC - Extract the PC context from a time-step context for PVode.
790: Input Parameter:
791: . ts - the time-step context
793: Output Parameter:
794: . pc - the preconditioner context
796: Level: advanced
798: Contributed by: Liyang Xu
800: .seealso: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
801: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
802: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
803: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance()
804: @*/
805: int TSPVodeGetPC(TS ts,PC *pc)
806: {
807: int ierr,(*f)(TS,PC *);
810: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeGetPC_C",(void (**)(void))&f);
811: if (f) {
812: (*f)(ts,pc);
813: } else {
814: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"TS must be of PVode type to extract the PC");
815: }
817: return(0);
818: }
822: /*@
823: TSPVodeSetExactFinalTime - Determines if PVode interpolates solution to the
824: exact final time requested by the user or just returns it at the final time
825: it computed. (Defaults to true).
827: Input Parameter:
828: + ts - the time-step context
829: - ft - PETSC_TRUE if interpolates, else PETSC_FALSE
831: Level: beginner
833: .seealso:TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
834: TSPVodeSetLinearTolerance(), TSPVodeSetGramSchmidtType(), TSPVodeSetTolerance(),
835: TSPVodeGetIterations(), TSPVodeSetType(), TSPVodeSetGMRESRestart(),
836: TSPVodeSetLinearTolerance(), TSPVodeSetTolerance(), TSPVodeGetPC()
837: @*/
838: int TSPVodeSetExactFinalTime(TS ts,PetscTruth ft)
839: {
840: int ierr,(*f)(TS,PetscTruth);
843: PetscObjectQueryFunction((PetscObject)ts,"TSPVodeSetExactFinalTime_C",(void (**)(void))&f);
844: if (f) {
845: (*f)(ts,ft);
846: }
848: return(0);
849: }
851: /* -------------------------------------------------------------------------------------------*/
853: /*
855: Contributed by: Liyang Xu
856: */
857: EXTERN_C_BEGIN
860: int TSCreate_PVode(TS ts)
861: {
862: TS_PVode *cvode;
863: int ierr;
866: ts->ops->destroy = TSDestroy_PVode;
867: ts->ops->view = TSView_PVode;
869: if (ts->problem_type != TS_NONLINEAR) {
870: SETERRQ(PETSC_ERR_SUP,"Only support for nonlinear problems");
871: }
872: ts->ops->setup = TSSetUp_PVode_Nonlinear;
873: ts->ops->step = TSStep_PVode_Nonlinear;
874: ts->ops->setfromoptions = TSSetFromOptions_PVode_Nonlinear;
876: PetscNew(TS_PVode,&cvode);
877: PetscMemzero(cvode,sizeof(TS_PVode));
878: PCCreate(ts->comm,&cvode->pc);
879: PetscLogObjectParent(ts,cvode->pc);
880: ts->data = (void*)cvode;
881: cvode->cvode_type = BDF;
882: cvode->gtype = PVODE_UNMODIFIED_GS;
883: cvode->restart = 5;
884: cvode->linear_tol = .05;
886: cvode->exact_final_time = PETSC_TRUE;
888: MPI_Comm_dup(ts->comm,&(cvode->comm_pvode));
889: /* set tolerance for PVode */
890: cvode->abstol = 1e-6;
891: cvode->reltol = 1e-6;
893: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetType_C","TSPVodeSetType_PVode",
894: TSPVodeSetType_PVode);
895: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGMRESRestart_C",
896: "TSPVodeSetGMRESRestart_PVode",
897: TSPVodeSetGMRESRestart_PVode);
898: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetLinearTolerance_C",
899: "TSPVodeSetLinearTolerance_PVode",
900: TSPVodeSetLinearTolerance_PVode);
901: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetGramSchmidtType_C",
902: "TSPVodeSetGramSchmidtType_PVode",
903: TSPVodeSetGramSchmidtType_PVode);
904: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetTolerance_C",
905: "TSPVodeSetTolerance_PVode",
906: TSPVodeSetTolerance_PVode);
907: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetPC_C",
908: "TSPVodeGetPC_PVode",
909: TSPVodeGetPC_PVode);
910: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeGetIterations_C",
911: "TSPVodeGetIterations_PVode",
912: TSPVodeGetIterations_PVode);
913: PetscObjectComposeFunctionDynamic((PetscObject)ts,"TSPVodeSetExactFinalTime_C",
914: "TSPVodeSetExactFinalTime_PVode",
915: TSPVodeSetExactFinalTime_PVode);
916: return(0);
917: }
918: EXTERN_C_END