Lines Matching refs:ts
67 PetscErrorCode TSIRKTableauCreate(TS ts, PetscInt nstages, const PetscReal *A, const PetscReal *b, … in TSIRKTableauCreate() argument
69 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKTableauCreate()
88 static PetscErrorCode TSIRKCreate_Gauss(TS ts) in TSIRKCreate_Gauss() argument
98 PetscCall(TSIRKGetNumStages(ts, &nstages)); in TSIRKCreate_Gauss()
153 …PetscCall(TSIRKTableauCreate(ts, nstages, gauss_A_real, gauss_b, gauss_c, NULL, gauss_A_inv, gauss… in TSIRKCreate_Gauss()
269 static PetscErrorCode TSEvaluateStep_IRK(TS ts, PetscInt order, Vec U, PetscBool *done) in TSEvaluateStep_IRK() argument
271 TS_IRK *irk = (TS_IRK *)ts->data; in TSEvaluateStep_IRK()
282 h = ts->time_step; in TSEvaluateStep_IRK()
285 h = ts->ptime - ts->ptime_prev; in TSEvaluateStep_IRK()
288 SETERRQ(PetscObjectComm((PetscObject)ts), PETSC_ERR_PLIB, "Invalid TSStepStatus"); in TSEvaluateStep_IRK()
291 PetscCall(VecCopy(ts->vec_sol, U)); in TSEvaluateStep_IRK()
297 static PetscErrorCode TSRollBack_IRK(TS ts) in TSRollBack_IRK() argument
299 TS_IRK *irk = (TS_IRK *)ts->data; in TSRollBack_IRK()
302 PetscCall(VecCopy(irk->U0, ts->vec_sol)); in TSRollBack_IRK()
306 static PetscErrorCode TSStep_IRK(TS ts) in TSStep_IRK() argument
308 TS_IRK *irk = (TS_IRK *)ts->data; in TSStep_IRK()
317 PetscReal next_time_step = ts->time_step; in TSStep_IRK()
320 if (!ts->steprollback) PetscCall(VecCopy(ts->vec_sol, irk->U0)); in TSStep_IRK()
321 PetscCall(VecGetBlockSize(ts->vec_sol, &bs)); in TSStep_IRK()
322 …for (i = 0; i < nstages; i++) PetscCall(VecStrideScatter(ts->vec_sol, i * bs, irk->Z, INSERT_VALUE… in TSStep_IRK()
325 while (!ts->reason && irk->status != TS_STEP_COMPLETE) { in TSStep_IRK()
326 PetscCall(VecCopy(ts->vec_sol, irk->U)); in TSStep_IRK()
327 PetscCall(TSGetSNES(ts, &snes)); in TSStep_IRK()
331 ts->snes_its += its; in TSStep_IRK()
332 ts->ksp_its += lits; in TSStep_IRK()
336 …for (j = 0; j < nstages; j++) PetscCall(VecAXPY(irk->YdotI[i], A_inv[i + j * nstages] / ts->time_s… in TSStep_IRK()
337 PetscCall(VecAXPY(irk->YdotI[i], -A_inv_rowsum[i] / ts->time_step, irk->U)); in TSStep_IRK()
340 PetscCall(TSEvaluateStep_IRK(ts, irk->order, ts->vec_sol, NULL)); in TSStep_IRK()
342 PetscCall(TSGetAdapt(ts, &adapt)); in TSStep_IRK()
343 PetscCall(TSAdaptChoose(adapt, ts, ts->time_step, NULL, &next_time_step, &accept)); in TSStep_IRK()
346 PetscCall(TSRollBack_IRK(ts)); in TSStep_IRK()
347 ts->time_step = next_time_step; in TSStep_IRK()
351 ts->ptime += ts->time_step; in TSStep_IRK()
352 ts->time_step = next_time_step; in TSStep_IRK()
355 ts->reject++; in TSStep_IRK()
357 if (!ts->reason && ++rejections > ts->max_reject && ts->max_reject >= 0) { in TSStep_IRK()
358 ts->reason = TS_DIVERGED_STEP_REJECTED; in TSStep_IRK()
359 …tscCall(PetscInfo(ts, "Step=%" PetscInt_FMT ", step rejections %" PetscInt_FMT " greater than curr… in TSStep_IRK()
365 static PetscErrorCode TSInterpolate_IRK(TS ts, PetscReal itime, Vec U) in TSInterpolate_IRK() argument
367 TS_IRK *irk = (TS_IRK *)ts->data; in TSInterpolate_IRK()
375 …PetscCheck(B, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "TSIRK %s does not have an interpol… in TSInterpolate_IRK()
379 h = ts->time_step; in TSInterpolate_IRK()
380 t = (itime - ts->ptime) / h; in TSInterpolate_IRK()
383 h = ts->ptime - ts->ptime_prev; in TSInterpolate_IRK()
384 t = (itime - ts->ptime) / h + 1; /* In the interval [0,1] */ in TSInterpolate_IRK()
387 SETERRQ(PetscObjectComm((PetscObject)ts), PETSC_ERR_PLIB, "Invalid TSStepStatus"); in TSInterpolate_IRK()
398 static PetscErrorCode TSIRKTableauReset(TS ts) in TSIRKTableauReset() argument
400 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKTableauReset()
410 static PetscErrorCode TSReset_IRK(TS ts) in TSReset_IRK() argument
412 TS_IRK *irk = (TS_IRK *)ts->data; in TSReset_IRK()
415 PetscCall(TSIRKTableauReset(ts)); in TSReset_IRK()
429 static PetscErrorCode TSIRKGetVecs(TS ts, DM dm, Vec *U) in TSIRKGetVecs() argument
431 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKGetVecs()
435 if (dm && dm != ts->dm) PetscCall(DMGetNamedGlobalVector(dm, "TSIRK_U", U)); in TSIRKGetVecs()
441 static PetscErrorCode TSIRKRestoreVecs(TS ts, DM dm, Vec *U) in TSIRKRestoreVecs() argument
445 if (dm && dm != ts->dm) PetscCall(DMRestoreNamedGlobalVector(dm, "TSIRK_U", U)); in TSIRKRestoreVecs()
456 static PetscErrorCode SNESTSFormFunction_IRK(SNES snes, Vec ZC, Vec FC, TS ts) in SNESTSFormFunction_IRK() argument
458 TS_IRK *irk = (TS_IRK *)ts->data; in SNESTSFormFunction_IRK()
465 PetscReal h = ts->time_step; in SNESTSFormFunction_IRK()
470 PetscCall(TSIRKGetVecs(ts, dm, &U)); in SNESTSFormFunction_IRK()
472 dmsave = ts->dm; in SNESTSFormFunction_IRK()
473 ts->dm = dm; in SNESTSFormFunction_IRK()
478 …PetscCall(TSComputeIFunction(ts, ts->ptime + ts->time_step * c[i], Y[i], Ydot, YdotI[i], PETSC_FAL… in SNESTSFormFunction_IRK()
481 ts->dm = dmsave; in SNESTSFormFunction_IRK()
482 PetscCall(TSIRKRestoreVecs(ts, dm, &U)); in SNESTSFormFunction_IRK()
492 static PetscErrorCode SNESTSFormJacobian_IRK(SNES snes, Vec ZC, Mat JC, Mat JCpre, TS ts) in SNESTSFormJacobian_IRK() argument
494 TS_IRK *irk = (TS_IRK *)ts->data; in SNESTSFormJacobian_IRK()
507 dmsave = ts->dm; in SNESTSFormJacobian_IRK()
508 ts->dm = dm; in SNESTSFormJacobian_IRK()
510 …PetscCheck(ts->equation_type <= TS_EQ_ODE_EXPLICIT, PetscObjectComm((PetscObject)ts), PETSC_ERR_SU… in SNESTSFormJacobian_IRK()
514 …PetscCall(TSComputeIJacobian(ts, ts->ptime + ts->time_step * c[nstages - 1], Y[nstages - 1], Ydot,… in SNESTSFormJacobian_IRK()
517 for (j = 0; j < nstages; j++) S[i + nstages * j] = tab->A_inv[i + nstages * j] / ts->time_step; in SNESTSFormJacobian_IRK()
519 ts->dm = dmsave; in SNESTSFormJacobian_IRK()
531 TS ts = (TS)ctx; in DMRestrictHook_TSIRK() local
535 PetscCall(TSIRKGetVecs(ts, fine, &U)); in DMRestrictHook_TSIRK()
536 PetscCall(TSIRKGetVecs(ts, coarse, &U_c)); in DMRestrictHook_TSIRK()
539 PetscCall(TSIRKRestoreVecs(ts, fine, &U)); in DMRestrictHook_TSIRK()
540 PetscCall(TSIRKRestoreVecs(ts, coarse, &U_c)); in DMRestrictHook_TSIRK()
552 TS ts = (TS)ctx; in DMSubDomainRestrictHook_TSIRK() local
556 PetscCall(TSIRKGetVecs(ts, dm, &U)); in DMSubDomainRestrictHook_TSIRK()
557 PetscCall(TSIRKGetVecs(ts, subdm, &U_c)); in DMSubDomainRestrictHook_TSIRK()
562 PetscCall(TSIRKRestoreVecs(ts, dm, &U)); in DMSubDomainRestrictHook_TSIRK()
563 PetscCall(TSIRKRestoreVecs(ts, subdm, &U_c)); in DMSubDomainRestrictHook_TSIRK()
567 static PetscErrorCode TSSetUp_IRK(TS ts) in TSSetUp_IRK() argument
569 TS_IRK *irk = (TS_IRK *)ts->data; in TSSetUp_IRK()
579 if (!irk->Y) PetscCall(VecDuplicateVecs(ts->vec_sol, irk->nstages, &irk->Y)); in TSSetUp_IRK()
580 if (!irk->YdotI) PetscCall(VecDuplicateVecs(ts->vec_sol, irk->nstages, &irk->YdotI)); in TSSetUp_IRK()
581 if (!irk->Ydot) PetscCall(VecDuplicate(ts->vec_sol, &irk->Ydot)); in TSSetUp_IRK()
582 if (!irk->U) PetscCall(VecDuplicate(ts->vec_sol, &irk->U)); in TSSetUp_IRK()
583 if (!irk->U0) PetscCall(VecDuplicate(ts->vec_sol, &irk->U0)); in TSSetUp_IRK()
585 PetscCall(VecCreate(PetscObjectComm((PetscObject)ts->vec_sol), &irk->Z)); in TSSetUp_IRK()
586 PetscCall(VecGetSize(ts->vec_sol, &vsize)); in TSSetUp_IRK()
588 PetscCall(VecGetBlockSize(ts->vec_sol, &bs)); in TSSetUp_IRK()
592 PetscCall(TSGetDM(ts, &dm)); in TSSetUp_IRK()
593 PetscCall(DMCoarsenHookAdd(dm, DMCoarsenHook_TSIRK, DMRestrictHook_TSIRK, ts)); in TSSetUp_IRK()
594 PetscCall(DMSubDomainHookAdd(dm, DMSubDomainHook_TSIRK, DMSubDomainRestrictHook_TSIRK, ts)); in TSSetUp_IRK()
596 PetscCall(TSGetSNES(ts, &ts->snes)); in TSSetUp_IRK()
598 PetscCall(SNESSetFunction(ts->snes, R, SNESTSFormFunction, ts)); in TSSetUp_IRK()
599 PetscCall(TSGetIJacobian(ts, &J, NULL, NULL, NULL)); in TSSetUp_IRK()
604 PetscCall(SNESSetJacobian(ts->snes, irk->TJ, irk->TJ, SNESTSFormJacobian, ts)); in TSSetUp_IRK()
609 static PetscErrorCode TSSetFromOptions_IRK(TS ts, PetscOptionItems PetscOptionsObject) in TSSetFromOptions_IRK() argument
611 TS_IRK *irk = (TS_IRK *)ts->data; in TSSetFromOptions_IRK()
621 PetscCall(TSIRKSetType(ts, tname)); in TSSetFromOptions_IRK()
628 static PetscErrorCode TSView_IRK(TS ts, PetscViewer viewer) in TSView_IRK() argument
630 TS_IRK *irk = (TS_IRK *)ts->data; in TSView_IRK()
640 PetscCall(TSIRKGetType(ts, &irktype)); in TSView_IRK()
651 static PetscErrorCode TSLoad_IRK(TS ts, PetscViewer viewer) in TSLoad_IRK() argument
657 PetscCall(TSGetAdapt(ts, &adapt)); in TSLoad_IRK()
659 PetscCall(TSGetSNES(ts, &snes)); in TSLoad_IRK()
662 PetscCall(SNESSetFunction(snes, NULL, NULL, ts)); in TSLoad_IRK()
663 PetscCall(SNESSetJacobian(snes, NULL, NULL, NULL, ts)); in TSLoad_IRK()
683 PetscErrorCode TSIRKSetType(TS ts, TSIRKType irktype) in TSIRKSetType() argument
686 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); in TSIRKSetType()
688 PetscTryMethod(ts, "TSIRKSetType_C", (TS, TSIRKType), (ts, irktype)); in TSIRKSetType()
707 PetscErrorCode TSIRKGetType(TS ts, TSIRKType *irktype) in TSIRKGetType() argument
710 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); in TSIRKGetType()
711 PetscUseMethod(ts, "TSIRKGetType_C", (TS, TSIRKType *), (ts, irktype)); in TSIRKGetType()
731 PetscErrorCode TSIRKSetNumStages(TS ts, PetscInt nstages) in TSIRKSetNumStages() argument
734 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); in TSIRKSetNumStages()
735 PetscTryMethod(ts, "TSIRKSetNumStages_C", (TS, PetscInt), (ts, nstages)); in TSIRKSetNumStages()
752 PetscErrorCode TSIRKGetNumStages(TS ts, PetscInt *nstages) in TSIRKGetNumStages() argument
755 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); in TSIRKGetNumStages()
757 PetscTryMethod(ts, "TSIRKGetNumStages_C", (TS, PetscInt *), (ts, nstages)); in TSIRKGetNumStages()
761 static PetscErrorCode TSIRKGetType_IRK(TS ts, TSIRKType *irktype) in TSIRKGetType_IRK() argument
763 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKGetType_IRK()
770 static PetscErrorCode TSIRKSetType_IRK(TS ts, TSIRKType irktype) in TSIRKSetType_IRK() argument
772 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKSetType_IRK()
778 PetscCall(TSIRKTableauReset(ts)); in TSIRKSetType_IRK()
781 …PetscCheck(irkcreate, PetscObjectComm((PetscObject)ts), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown TSIRK… in TSIRKSetType_IRK()
782 PetscCall((*irkcreate)(ts)); in TSIRKSetType_IRK()
787 static PetscErrorCode TSIRKSetNumStages_IRK(TS ts, PetscInt nstages) in TSIRKSetNumStages_IRK() argument
789 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKSetNumStages_IRK()
797 static PetscErrorCode TSIRKGetNumStages_IRK(TS ts, PetscInt *nstages) in TSIRKGetNumStages_IRK() argument
799 TS_IRK *irk = (TS_IRK *)ts->data; in TSIRKGetNumStages_IRK()
807 static PetscErrorCode TSDestroy_IRK(TS ts) in TSDestroy_IRK() argument
810 PetscCall(TSReset_IRK(ts)); in TSDestroy_IRK()
811 if (ts->dm) { in TSDestroy_IRK()
812 PetscCall(DMCoarsenHookRemove(ts->dm, DMCoarsenHook_TSIRK, DMRestrictHook_TSIRK, ts)); in TSDestroy_IRK()
813 …PetscCall(DMSubDomainHookRemove(ts->dm, DMSubDomainHook_TSIRK, DMSubDomainRestrictHook_TSIRK, ts)); in TSDestroy_IRK()
815 PetscCall(PetscFree(ts->data)); in TSDestroy_IRK()
816 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKSetType_C", NULL)); in TSDestroy_IRK()
817 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKGetType_C", NULL)); in TSDestroy_IRK()
818 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKSetNumStages_C", NULL)); in TSDestroy_IRK()
819 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKGetNumStages_C", NULL)); in TSDestroy_IRK()
837 PETSC_EXTERN PetscErrorCode TSCreate_IRK(TS ts) in TSCreate_IRK() argument
844 ts->ops->reset = TSReset_IRK; in TSCreate_IRK()
845 ts->ops->destroy = TSDestroy_IRK; in TSCreate_IRK()
846 ts->ops->view = TSView_IRK; in TSCreate_IRK()
847 ts->ops->load = TSLoad_IRK; in TSCreate_IRK()
848 ts->ops->setup = TSSetUp_IRK; in TSCreate_IRK()
849 ts->ops->step = TSStep_IRK; in TSCreate_IRK()
850 ts->ops->interpolate = TSInterpolate_IRK; in TSCreate_IRK()
851 ts->ops->evaluatestep = TSEvaluateStep_IRK; in TSCreate_IRK()
852 ts->ops->rollback = TSRollBack_IRK; in TSCreate_IRK()
853 ts->ops->setfromoptions = TSSetFromOptions_IRK; in TSCreate_IRK()
854 ts->ops->snesfunction = SNESTSFormFunction_IRK; in TSCreate_IRK()
855 ts->ops->snesjacobian = SNESTSFormJacobian_IRK; in TSCreate_IRK()
857 ts->usessnes = PETSC_TRUE; in TSCreate_IRK()
860 ts->data = (void *)irk; in TSCreate_IRK()
862 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKSetType_C", TSIRKSetType_IRK)); in TSCreate_IRK()
863 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKGetType_C", TSIRKGetType_IRK)); in TSCreate_IRK()
864 …PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKSetNumStages_C", TSIRKSetNumStages_IRK… in TSCreate_IRK()
865 …PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSIRKGetNumStages_C", TSIRKGetNumStages_IRK… in TSCreate_IRK()
869 PetscCall(TSIRKSetType(ts, TSIRKDefault)); in TSCreate_IRK()