1 /* 2 Code for timestepping with discrete gradient integrators 3 */ 4 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 5 #include <petscdm.h> 6 7 PetscBool DGCite = PETSC_FALSE; 8 const char DGCitation[] = "@article{Gonzalez1996,\n" 9 " title = {Time integration and discrete Hamiltonian systems},\n" 10 " author = {Oscar Gonzalez},\n" 11 " journal = {Journal of Nonlinear Science},\n" 12 " volume = {6},\n" 13 " pages = {449--467},\n" 14 " doi = {10.1007/978-1-4612-1246-1_10},\n" 15 " year = {1996}\n}\n"; 16 17 typedef struct { 18 PetscReal stage_time; 19 Vec X0, X, Xdot; 20 PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *); 21 PetscErrorCode (*Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *); 22 PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *); 23 } TS_DiscGrad; 24 25 static PetscErrorCode TSDiscGradGetX0AndXdot(TS ts, DM dm, Vec *X0, Vec *Xdot) 26 { 27 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 28 PetscErrorCode ierr; 29 30 PetscFunctionBegin; 31 if (X0) { 32 if (dm && dm != ts->dm) {ierr = DMGetNamedGlobalVector(dm, "TSDiscGrad_X0", X0);CHKERRQ(ierr);} 33 else {*X0 = ts->vec_sol;} 34 } 35 if (Xdot) { 36 if (dm && dm != ts->dm) {ierr = DMGetNamedGlobalVector(dm, "TSDiscGrad_Xdot", Xdot);CHKERRQ(ierr);} 37 else {*Xdot = dg->Xdot;} 38 } 39 PetscFunctionReturn(0); 40 } 41 42 static PetscErrorCode TSDiscGradRestoreX0AndXdot(TS ts, DM dm, Vec *X0, Vec *Xdot) 43 { 44 PetscErrorCode ierr; 45 46 PetscFunctionBegin; 47 if (X0) { 48 if (dm && dm != ts->dm) {ierr = DMRestoreNamedGlobalVector(dm, "TSDiscGrad_X0", X0);CHKERRQ(ierr);} 49 } 50 if (Xdot) { 51 if (dm && dm != ts->dm) {ierr = DMRestoreNamedGlobalVector(dm, "TSDiscGrad_Xdot", Xdot);CHKERRQ(ierr);} 52 } 53 PetscFunctionReturn(0); 54 } 55 56 static PetscErrorCode DMCoarsenHook_TSDiscGrad(DM fine, DM coarse, void *ctx) 57 { 58 PetscFunctionBegin; 59 PetscFunctionReturn(0); 60 } 61 62 static PetscErrorCode DMRestrictHook_TSDiscGrad(DM fine, Mat restrct, Vec rscale, Mat inject, DM coarse, void *ctx) 63 { 64 TS ts = (TS) ctx; 65 Vec X0, Xdot, X0_c, Xdot_c; 66 PetscErrorCode ierr; 67 68 PetscFunctionBegin; 69 ierr = TSDiscGradGetX0AndXdot(ts, fine, &X0, &Xdot);CHKERRQ(ierr); 70 ierr = TSDiscGradGetX0AndXdot(ts, coarse, &X0_c, &Xdot_c);CHKERRQ(ierr); 71 ierr = MatRestrict(restrct, X0, X0_c);CHKERRQ(ierr); 72 ierr = MatRestrict(restrct, Xdot, Xdot_c);CHKERRQ(ierr); 73 ierr = VecPointwiseMult(X0_c, rscale, X0_c);CHKERRQ(ierr); 74 ierr = VecPointwiseMult(Xdot_c, rscale, Xdot_c);CHKERRQ(ierr); 75 ierr = TSDiscGradRestoreX0AndXdot(ts, fine, &X0, &Xdot);CHKERRQ(ierr); 76 ierr = TSDiscGradRestoreX0AndXdot(ts, coarse, &X0_c, &Xdot_c);CHKERRQ(ierr); 77 PetscFunctionReturn(0); 78 } 79 80 static PetscErrorCode DMSubDomainHook_TSDiscGrad(DM dm, DM subdm, void *ctx) 81 { 82 PetscFunctionBegin; 83 PetscFunctionReturn(0); 84 } 85 86 static PetscErrorCode DMSubDomainRestrictHook_TSDiscGrad(DM dm, VecScatter gscat, VecScatter lscat, DM subdm, void *ctx) 87 { 88 TS ts = (TS) ctx; 89 Vec X0, Xdot, X0_sub, Xdot_sub; 90 PetscErrorCode ierr; 91 92 PetscFunctionBegin; 93 ierr = TSDiscGradGetX0AndXdot(ts, dm, &X0, &Xdot);CHKERRQ(ierr); 94 ierr = TSDiscGradGetX0AndXdot(ts, subdm, &X0_sub, &Xdot_sub);CHKERRQ(ierr); 95 96 ierr = VecScatterBegin(gscat, X0, X0_sub, INSERT_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 97 ierr = VecScatterEnd(gscat, X0, X0_sub, INSERT_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 98 99 ierr = VecScatterBegin(gscat, Xdot, Xdot_sub, INSERT_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 100 ierr = VecScatterEnd(gscat, Xdot, Xdot_sub, INSERT_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 101 102 ierr = TSDiscGradRestoreX0AndXdot(ts, dm, &X0, &Xdot);CHKERRQ(ierr); 103 ierr = TSDiscGradRestoreX0AndXdot(ts, subdm, &X0_sub, &Xdot_sub);CHKERRQ(ierr); 104 PetscFunctionReturn(0); 105 } 106 107 static PetscErrorCode TSSetUp_DiscGrad(TS ts) 108 { 109 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 110 DM dm; 111 PetscErrorCode ierr; 112 113 PetscFunctionBegin; 114 if (!dg->X) {ierr = VecDuplicate(ts->vec_sol, &dg->X);CHKERRQ(ierr);} 115 if (!dg->X0) {ierr = VecDuplicate(ts->vec_sol, &dg->X0);CHKERRQ(ierr);} 116 if (!dg->Xdot) {ierr = VecDuplicate(ts->vec_sol, &dg->Xdot);CHKERRQ(ierr);} 117 118 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 119 ierr = DMCoarsenHookAdd(dm, DMCoarsenHook_TSDiscGrad, DMRestrictHook_TSDiscGrad, ts);CHKERRQ(ierr); 120 ierr = DMSubDomainHookAdd(dm, DMSubDomainHook_TSDiscGrad, DMSubDomainRestrictHook_TSDiscGrad, ts);CHKERRQ(ierr); 121 PetscFunctionReturn(0); 122 } 123 124 static PetscErrorCode TSSetFromOptions_DiscGrad(PetscOptionItems *PetscOptionsObject, TS ts) 125 { 126 PetscErrorCode ierr; 127 128 PetscFunctionBegin; 129 ierr = PetscOptionsHead(PetscOptionsObject, "Discrete Gradients ODE solver options");CHKERRQ(ierr); 130 { 131 //ierr = PetscOptionsReal("-ts_theta_theta","Location of stage (0<Theta<=1)","TSDiscGradSetTheta",th->Theta,&th->Theta,NULL);CHKERRQ(ierr); 132 } 133 ierr = PetscOptionsTail();CHKERRQ(ierr); 134 PetscFunctionReturn(0); 135 } 136 137 static PetscErrorCode TSView_DiscGrad(TS ts,PetscViewer viewer) 138 { 139 PetscBool iascii; 140 PetscErrorCode ierr; 141 142 PetscFunctionBegin; 143 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 144 if (iascii) { 145 ierr = PetscViewerASCIIPrintf(viewer," Discrete Gradients\n");CHKERRQ(ierr); 146 } 147 PetscFunctionReturn(0); 148 } 149 150 static PetscErrorCode TSReset_DiscGrad(TS ts) 151 { 152 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 153 PetscErrorCode ierr; 154 155 PetscFunctionBegin; 156 ierr = VecDestroy(&dg->X);CHKERRQ(ierr); 157 ierr = VecDestroy(&dg->X0);CHKERRQ(ierr); 158 ierr = VecDestroy(&dg->Xdot);CHKERRQ(ierr); 159 PetscFunctionReturn(0); 160 } 161 162 static PetscErrorCode TSDestroy_DiscGrad(TS ts) 163 { 164 DM dm; 165 PetscErrorCode ierr; 166 167 PetscFunctionBegin; 168 ierr = TSReset_DiscGrad(ts);CHKERRQ(ierr); 169 ierr = TSGetDM(ts, &dm);CHKERRQ(ierr); 170 if (dm) { 171 ierr = DMCoarsenHookRemove(dm, DMCoarsenHook_TSDiscGrad, DMRestrictHook_TSDiscGrad, ts);CHKERRQ(ierr); 172 ierr = DMSubDomainHookRemove(dm, DMSubDomainHook_TSDiscGrad, DMSubDomainRestrictHook_TSDiscGrad, ts);CHKERRQ(ierr); 173 } 174 ierr = PetscFree(ts->data);CHKERRQ(ierr); 175 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradGetFormulation_C",NULL);CHKERRQ(ierr); 176 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradSetFormulation_C",NULL);CHKERRQ(ierr); 177 PetscFunctionReturn(0); 178 } 179 180 static PetscErrorCode TSInterpolate_DiscGrad(TS ts, PetscReal t, Vec X) 181 { 182 TS_DiscGrad *dg = (TS_DiscGrad*)ts->data; 183 PetscReal dt = t - ts->ptime; 184 PetscErrorCode ierr; 185 186 PetscFunctionBegin; 187 ierr = VecCopy(ts->vec_sol, dg->X);CHKERRQ(ierr); 188 ierr = VecWAXPY(X, dt, dg->Xdot, dg->X);CHKERRQ(ierr); 189 PetscFunctionReturn(0); 190 } 191 192 static PetscErrorCode TSDiscGrad_SNESSolve(TS ts, Vec b, Vec x) 193 { 194 SNES snes; 195 PetscInt nits, lits; 196 PetscErrorCode ierr; 197 198 PetscFunctionBegin; 199 ierr = TSGetSNES(ts, &snes);CHKERRQ(ierr); 200 ierr = SNESSolve(snes, b, x);CHKERRQ(ierr); 201 ierr = SNESGetIterationNumber(snes, &nits);CHKERRQ(ierr); 202 ierr = SNESGetLinearSolveIterations(snes, &lits);CHKERRQ(ierr); 203 ts->snes_its += nits; 204 ts->ksp_its += lits; 205 PetscFunctionReturn(0); 206 } 207 208 static PetscErrorCode TSStep_DiscGrad(TS ts) 209 { 210 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 211 TSAdapt adapt; 212 TSStepStatus status = TS_STEP_INCOMPLETE; 213 PetscInt rejections = 0; 214 PetscBool stageok, accept = PETSC_TRUE; 215 PetscReal next_time_step = ts->time_step; 216 PetscErrorCode ierr; 217 218 PetscFunctionBegin; 219 ierr = TSGetAdapt(ts, &adapt);CHKERRQ(ierr); 220 if (!ts->steprollback) {ierr = VecCopy(ts->vec_sol, dg->X0);CHKERRQ(ierr);} 221 222 while (!ts->reason && status != TS_STEP_COMPLETE) { 223 PetscReal shift = 1/(0.5*ts->time_step); 224 225 dg->stage_time = ts->ptime + 0.5*ts->time_step; 226 227 ierr = VecCopy(dg->X0, dg->X);CHKERRQ(ierr); 228 ierr = TSPreStage(ts, dg->stage_time);CHKERRQ(ierr); 229 ierr = TSDiscGrad_SNESSolve(ts, NULL, dg->X);CHKERRQ(ierr); 230 ierr = TSPostStage(ts, dg->stage_time, 0, &dg->X);CHKERRQ(ierr); 231 ierr = TSAdaptCheckStage(adapt, ts, dg->stage_time, dg->X, &stageok);CHKERRQ(ierr); 232 if (!stageok) goto reject_step; 233 234 status = TS_STEP_PENDING; 235 ierr = VecAXPBYPCZ(dg->Xdot, -shift, shift, 0, dg->X0, dg->X);CHKERRQ(ierr); 236 ierr = VecAXPY(ts->vec_sol, ts->time_step, dg->Xdot);CHKERRQ(ierr); 237 ierr = TSAdaptChoose(adapt, ts, ts->time_step, NULL, &next_time_step, &accept);CHKERRQ(ierr); 238 status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE; 239 if (!accept) { 240 ierr = VecCopy(dg->X0, ts->vec_sol);CHKERRQ(ierr); 241 ts->time_step = next_time_step; 242 goto reject_step; 243 } 244 ts->ptime += ts->time_step; 245 ts->time_step = next_time_step; 246 break; 247 248 reject_step: 249 ts->reject++; accept = PETSC_FALSE; 250 if (!ts->reason && ts->max_reject >= 0 && ++rejections > ts->max_reject) { 251 ts->reason = TS_DIVERGED_STEP_REJECTED; 252 ierr = PetscInfo2(ts, "Step=%D, step rejections %D greater than current TS allowed, stopping solve\n", ts->steps, rejections);CHKERRQ(ierr); 253 } 254 } 255 PetscFunctionReturn(0); 256 } 257 258 static PetscErrorCode TSGetStages_DiscGrad(TS ts, PetscInt *ns, Vec **Y) 259 { 260 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 261 262 PetscFunctionBegin; 263 if (ns) *ns = 1; 264 if (Y) *Y = &(dg->X); 265 PetscFunctionReturn(0); 266 } 267 268 /* 269 This defines the nonlinear equation that is to be solved with SNES 270 G(U) = F[t0 + 0.5*dt, U, (U-U0)/dt] = 0 271 */ 272 static PetscErrorCode SNESTSFormFunction_DiscGrad(SNES snes, Vec x, Vec y, TS ts) 273 { 274 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 275 PetscReal shift = 1/(0.5*ts->time_step); 276 Vec X0, Xdot; 277 DM dm, dmsave; 278 PetscErrorCode ierr; 279 280 PetscFunctionBegin; 281 ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr); 282 ierr = TSDiscGradGetX0AndXdot(ts, dm, &X0, &Xdot);CHKERRQ(ierr); 283 ierr = VecAXPBYPCZ(Xdot, -shift, shift, 0, X0, x);CHKERRQ(ierr); 284 285 /* DM monkey-business allows user code to call TSGetDM() inside of functions evaluated on levels of FAS */ 286 dmsave = ts->dm; 287 ts->dm = dm; 288 ierr = TSComputeIFunction(ts, dg->stage_time, x, Xdot, y, PETSC_FALSE);CHKERRQ(ierr); 289 ts->dm = dmsave; 290 ierr = TSDiscGradRestoreX0AndXdot(ts, dm, &X0, &Xdot);CHKERRQ(ierr); 291 PetscFunctionReturn(0); 292 } 293 294 static PetscErrorCode SNESTSFormJacobian_DiscGrad(SNES snes, Vec x, Mat A, Mat B, TS ts) 295 { 296 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 297 PetscReal shift = 1/(0.5*ts->time_step); 298 Vec Xdot; 299 DM dm,dmsave; 300 PetscErrorCode ierr; 301 302 PetscFunctionBegin; 303 ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr); 304 /* Xdot has already been computed in SNESTSFormFunction_DiscGrad (SNES guarantees this) */ 305 ierr = TSDiscGradGetX0AndXdot(ts, dm, NULL, &Xdot);CHKERRQ(ierr); 306 307 dmsave = ts->dm; 308 ts->dm = dm; 309 ierr = TSComputeIJacobian(ts, dg->stage_time, x, Xdot, shift, A, B, PETSC_FALSE);CHKERRQ(ierr); 310 ts->dm = dmsave; 311 ierr = TSDiscGradRestoreX0AndXdot(ts, dm, NULL, &Xdot);CHKERRQ(ierr); 312 PetscFunctionReturn(0); 313 } 314 315 static PetscErrorCode TSDiscGradGetFormulation_DiscGrad(TS ts, PetscErrorCode (**Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (**Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (**Gfunc)(TS, PetscReal, Vec, Vec, void *)) 316 { 317 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 318 319 PetscFunctionBegin; 320 *Sfunc = dg->Sfunc; 321 *Ffunc = dg->Ffunc; 322 *Gfunc = dg->Gfunc; 323 PetscFunctionReturn(0); 324 } 325 326 static PetscErrorCode TSDiscGradSetFormulation_DiscGrad(TS ts, PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (*Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *)) 327 { 328 TS_DiscGrad *dg = (TS_DiscGrad *) ts->data; 329 330 PetscFunctionBegin; 331 dg->Sfunc = Sfunc; 332 dg->Ffunc = Ffunc; 333 dg->Gfunc = Gfunc; 334 PetscFunctionReturn(0); 335 } 336 337 /*MC 338 TSDISCGRAD - ODE solver using the discrete gradients version of the implicit midpoint method 339 340 Level: beginner 341 342 Notes: This is the implicit midpoint rule, with an optional term that guarantees the discrete gradient property. This 343 timestepper applies to systems of the form 344 $ u_t = S(u) grad F(u) 345 where S(u) is a linear operator, and F is a functional of u. 346 347 .seealso: TSCreate(), TSSetType(), TS, TSTHETA, TSDiscGradSetFormulation() 348 M*/ 349 PETSC_EXTERN PetscErrorCode TSCreate_DiscGrad(TS ts) 350 { 351 TS_DiscGrad *th; 352 PetscErrorCode ierr; 353 354 PetscFunctionBegin; 355 ierr = PetscCitationsRegister(DGCitation, &DGCite);CHKERRQ(ierr); 356 ts->ops->reset = TSReset_DiscGrad; 357 ts->ops->destroy = TSDestroy_DiscGrad; 358 ts->ops->view = TSView_DiscGrad; 359 ts->ops->setfromoptions = TSSetFromOptions_DiscGrad; 360 ts->ops->setup = TSSetUp_DiscGrad; 361 ts->ops->step = TSStep_DiscGrad; 362 ts->ops->interpolate = TSInterpolate_DiscGrad; 363 ts->ops->getstages = TSGetStages_DiscGrad; 364 ts->ops->snesfunction = SNESTSFormFunction_DiscGrad; 365 ts->ops->snesjacobian = SNESTSFormJacobian_DiscGrad; 366 ts->default_adapt_type = TSADAPTNONE; 367 368 ts->usessnes = PETSC_TRUE; 369 370 ierr = PetscNewLog(ts,&th);CHKERRQ(ierr); 371 ts->data = (void*)th; 372 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradGetFormulation_C",TSDiscGradGetFormulation_DiscGrad);CHKERRQ(ierr); 373 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradSetFormulation_C",TSDiscGradSetFormulation_DiscGrad);CHKERRQ(ierr); 374 PetscFunctionReturn(0); 375 } 376 377 /*@C 378 TSDiscGradGetFormulation - Get the construction method for S, F, and grad F from the formulation u_t = S grad F 379 380 Not Collective 381 382 Input Parameter: 383 . ts - timestepping context 384 385 Output Parameters: 386 + Sfunc - constructor for the S matrix from the formulation 387 . Ffunc - functional F from the formulation 388 - Gfunc - constructor for the gradient of F from the formulation 389 390 Calling sequence of Sfunc: 391 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Mat S, void *) 392 393 Calling sequence of Ffunc: 394 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, PetscScalar *F, void *) 395 396 Calling sequence of Gfunc: 397 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Vec G, void *) 398 399 Level: intermediate 400 401 .seealso: TSDiscGradSetFormulation() 402 @*/ 403 PetscErrorCode TSDiscGradGetFormulation(TS ts, PetscErrorCode (**Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (**Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (**Gfunc)(TS, PetscReal, Vec, Vec, void *)) 404 { 405 PetscErrorCode ierr; 406 407 PetscFunctionBegin; 408 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 409 PetscValidPointer(Sfunc, 2); 410 PetscValidPointer(Ffunc, 3); 411 PetscValidPointer(Gfunc, 4); 412 ierr = PetscUseMethod(ts,"TSDiscGradGetFormulation_C",(TS,PetscErrorCode(**Sfunc)(TS,PetscReal,Vec,Mat,void*),PetscErrorCode(**Ffunc)(TS,PetscReal,Vec,PetscScalar*,void*),PetscErrorCode(**Gfunc)(TS,PetscReal,Vec,Vec,void*)),(ts,Sfunc,Ffunc,Gfunc));CHKERRQ(ierr); 413 PetscFunctionReturn(0); 414 } 415 416 /*@C 417 TSDiscGradSetFormulation - Set the construction method for S, F, and grad F from the formulation u_t = S(u) grad F(u) 418 419 Not Collective 420 421 Input Parameters: 422 + ts - timestepping context 423 . Sfunc - constructor for the S matrix from the formulation 424 . Ffunc - functional F from the formulation 425 - Gfunc - constructor for the gradient of F from the formulation 426 427 Calling sequence of Sfunc: 428 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Mat S, void *) 429 430 Calling sequence of Ffunc: 431 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, PetscScalar *F, void *) 432 433 Calling sequence of Gfunc: 434 $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Vec G, void *) 435 436 Level: Intermediate 437 438 .seealso: TSDiscGradGetFormulation() 439 @*/ 440 PetscErrorCode TSDiscGradSetFormulation(TS ts, PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (*Ffunc)(TS, PetscReal, Vec , PetscScalar *, void *), PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *), void *ctx) 441 { 442 PetscErrorCode ierr; 443 444 PetscFunctionBegin; 445 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 446 PetscValidFunction(Sfunc, 2); 447 PetscValidFunction(Ffunc, 3); 448 PetscValidFunction(Gfunc, 4); 449 ierr = PetscTryMethod(ts,"TSDiscGradSetFormulation_C",(TS,PetscErrorCode(*Sfunc)(TS,PetscReal,Vec,Mat,void*),PetscErrorCode(*Ffunc)(TS,PetscReal,Vec,PetscScalar*,void*),PetscErrorCode(*Gfunc)(TS,PetscReal,Vec,Vec,void*)),(ts,Sfunc,Ffunc,Gfunc));CHKERRQ(ierr); 450 PetscFunctionReturn(0); 451 } 452