Lines Matching refs:th
47 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaGetX0AndXdot() local
56 else *Xdot = th->Xdot; in TSThetaGetX0AndXdot()
124 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaEvaluateCostIntegral() local
128 if (th->endpoint) { in TSThetaEvaluateCostIntegral()
130 if (th->Theta != 1.0) { in TSThetaEvaluateCostIntegral()
131 PetscCall(TSComputeRHSFunction(quadts, th->ptime0, th->X0, ts->vec_costintegrand)); in TSThetaEvaluateCostIntegral()
132 … PetscCall(VecAXPY(quadts->vec_sol, th->time_step0 * (1.0 - th->Theta), ts->vec_costintegrand)); in TSThetaEvaluateCostIntegral()
135 PetscCall(VecAXPY(quadts->vec_sol, th->time_step0 * th->Theta, ts->vec_costintegrand)); in TSThetaEvaluateCostIntegral()
137 PetscCall(TSComputeRHSFunction(quadts, th->stage_time, th->X, ts->vec_costintegrand)); in TSThetaEvaluateCostIntegral()
138 PetscCall(VecAXPY(quadts->vec_sol, th->time_step0, ts->vec_costintegrand)); in TSThetaEvaluateCostIntegral()
145 TS_Theta *th = (TS_Theta *)ts->data; in TSForwardCostIntegral_Theta() local
150 PetscCall(VecCopy(quadts->vec_sol, th->VecCostIntegral0)); in TSForwardCostIntegral_Theta()
157 TS_Theta *th = (TS_Theta *)ts->data; in TSAdjointCostIntegral_Theta() local
161 th->ptime0 = ts->ptime + ts->time_step; in TSAdjointCostIntegral_Theta()
162 th->time_step0 = -ts->time_step; in TSAdjointCostIntegral_Theta()
182 TS_Theta *th = (TS_Theta *)ts->data; in TSResizeRegister_Theta() local
186 PetscCall(TSResizeRegisterVec(ts, "ts:theta:sol_prev", th->vec_sol_prev)); in TSResizeRegister_Theta()
187 PetscCall(TSResizeRegisterVec(ts, "ts:theta:X0", th->X0)); in TSResizeRegister_Theta()
189 PetscCall(TSResizeRetrieveVec(ts, "ts:theta:sol_prev", &th->vec_sol_prev)); in TSResizeRegister_Theta()
190 PetscCall(PetscObjectReference((PetscObject)th->vec_sol_prev)); in TSResizeRegister_Theta()
191 PetscCall(TSResizeRetrieveVec(ts, "ts:theta:X0", &th->X0)); in TSResizeRegister_Theta()
192 PetscCall(PetscObjectReference((PetscObject)th->X0)); in TSResizeRegister_Theta()
199 TS_Theta *th = (TS_Theta *)ts->data; in TSStep_Theta() local
206 if (th->vec_sol_prev) PetscCall(VecCopy(th->X0, th->vec_sol_prev)); in TSStep_Theta()
207 PetscCall(VecCopy(ts->vec_sol, th->X0)); in TSStep_Theta()
210 th->status = TS_STEP_INCOMPLETE; in TSStep_Theta()
211 while (!ts->reason && th->status != TS_STEP_COMPLETE) { in TSStep_Theta()
212 th->shift = 1 / (th->Theta * ts->time_step); in TSStep_Theta()
213 th->stage_time = ts->ptime + (th->endpoint ? (PetscReal)1 : th->Theta) * ts->time_step; in TSStep_Theta()
214 PetscCall(VecCopy(th->X0, th->X)); in TSStep_Theta()
215 if (th->extrapolate && !ts->steprestart) PetscCall(VecAXPY(th->X, 1 / th->shift, th->Xdot)); in TSStep_Theta()
216 if (th->endpoint) { /* This formulation assumes linear time-independent mass matrix */ in TSStep_Theta()
217 if (!th->affine) PetscCall(VecDuplicate(ts->vec_sol, &th->affine)); in TSStep_Theta()
218 PetscCall(VecZeroEntries(th->Xdot)); in TSStep_Theta()
219 PetscCall(TSComputeIFunction(ts, ts->ptime, th->X0, th->Xdot, th->affine, PETSC_FALSE)); in TSStep_Theta()
220 PetscCall(VecScale(th->affine, (th->Theta - 1) / th->Theta)); in TSStep_Theta()
222 PetscCall(TSPreStage(ts, th->stage_time)); in TSStep_Theta()
223 PetscCall(TSTheta_SNESSolve(ts, th->endpoint ? th->affine : NULL, th->X)); in TSStep_Theta()
224 PetscCall(TSPostStage(ts, th->stage_time, 0, &th->X)); in TSStep_Theta()
225 PetscCall(TSAdaptCheckStage(ts->adapt, ts, th->stage_time, th->X, &stageok)); in TSStep_Theta()
228 if (th->endpoint) { in TSStep_Theta()
229 PetscCall(VecCopy(th->X, ts->vec_sol)); in TSStep_Theta()
231 …PetscCall(VecAXPBYPCZ(th->Xdot, -th->shift, th->shift, 0, th->X0, th->X)); /* th->Xdot is needed b… in TSStep_Theta()
232 …if (th->Theta == 1.0) PetscCall(VecCopy(th->X, ts->vec_sol)); /* BEULER, stage alread… in TSStep_Theta()
234 PetscCall(VecAXPY(ts->vec_sol, ts->time_step, th->Xdot)); in TSStep_Theta()
237 PetscCall(VecCopy(th->X0, ts->vec_sol)); in TSStep_Theta()
243 th->status = TS_STEP_PENDING; in TSStep_Theta()
245 th->status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE; in TSStep_Theta()
247 PetscCall(VecCopy(th->X0, ts->vec_sol)); in TSStep_Theta()
253 th->ptime0 = ts->ptime; in TSStep_Theta()
254 th->time_step0 = ts->time_step; in TSStep_Theta()
273 TS_Theta *th = (TS_Theta *)ts->data; in TSAdjointStepBEuler_Private() local
275 …Vec *VecsDeltaLam = th->VecsDeltaLam, *VecsDeltaMu = th->VecsDeltaMu, *VecsSensiTemp = t… in TSAdjointStepBEuler_Private()
276 …Vec *VecsDeltaLam2 = th->VecsDeltaLam2, *VecsDeltaMu2 = th->VecsDeltaMu2, *VecsSensi2Tem… in TSAdjointStepBEuler_Private()
292 th->status = TS_STEP_INCOMPLETE; in TSAdjointStepBEuler_Private()
300 th->stage_time = ts->ptime; in TSAdjointStepBEuler_Private()
304 if (quadts) PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, ts->vec_sol, quadJ, NULL)); in TSAdjointStepBEuler_Private()
319 th->shift = 1. / adjoint_time_step; in TSAdjointStepBEuler_Private()
339 …PetscCall(TSComputeIHessianProductFunctionUU(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStepBEuler_Private()
341 …PetscCall(TSComputeIHessianProductFunctionUP(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStepBEuler_Private()
362 th->shift = 0.0; in TSAdjointStepBEuler_Private()
379 …PetscCall(VecAXPBYPCZ(th->Xdot, -1. / adjoint_time_step, 1.0 / adjoint_time_step, 0, th->X0, ts->v… in TSAdjointStepBEuler_Private()
380 …PetscCall(TSComputeIJacobianP(ts, th->stage_time, ts->vec_sol, th->Xdot, 1. / adjoint_time_step, t… in TSAdjointStepBEuler_Private()
381 if (quadts) PetscCall(TSComputeRHSJacobianP(quadts, th->stage_time, ts->vec_sol, quadJp)); in TSAdjointStepBEuler_Private()
384 …PetscCall(TSComputeIHessianProductFunctionPU(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStepBEuler_Private()
386 …PetscCall(TSComputeIHessianProductFunctionPP(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStepBEuler_Private()
412 th->status = TS_STEP_COMPLETE; in TSAdjointStepBEuler_Private()
418 TS_Theta *th = (TS_Theta *)ts->data; in TSAdjointStep_Theta() local
420 …Vec *VecsDeltaLam = th->VecsDeltaLam, *VecsDeltaMu = th->VecsDeltaMu, *VecsSensiTemp = th-… in TSAdjointStep_Theta()
421 …Vec *VecsDeltaLam2 = th->VecsDeltaLam2, *VecsDeltaMu2 = th->VecsDeltaMu2, *VecsSensi2Temp … in TSAdjointStep_Theta()
430 if (th->Theta == 1.) { in TSAdjointStep_Theta()
434 th->status = TS_STEP_INCOMPLETE; in TSAdjointStep_Theta()
442 th->stage_time = th->endpoint ? ts->ptime : (ts->ptime + (1. - th->Theta) * ts->time_step); in TSAdjointStep_Theta()
446 if (!th->endpoint) { in TSAdjointStep_Theta()
448 …PetscCall(VecAXPBYPCZ(th->X0, 1.0 / (1.0 - th->Theta), th->Theta / (th->Theta - 1.0), 0, th->X, ts… in TSAdjointStep_Theta()
454 if (th->endpoint) { in TSAdjointStep_Theta()
455 PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, ts->vec_sol, quadJ, NULL)); in TSAdjointStep_Theta()
457 PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, th->X, quadJ, NULL)); in TSAdjointStep_Theta()
463 PetscCall(VecScale(VecsSensiTemp[nadj], 1. / (th->Theta * adjoint_time_step))); in TSAdjointStep_Theta()
474 th->shift = 1. / (th->Theta * adjoint_time_step); in TSAdjointStep_Theta()
475 if (th->endpoint) { in TSAdjointStep_Theta()
478 PetscCall(TSComputeSNESJacobian(ts, th->X, J, Jpre)); in TSAdjointStep_Theta()
495 …PetscCheck(th->endpoint, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Operation not implement… in TSAdjointStep_Theta()
500 …PetscCall(TSComputeIHessianProductFunctionUU(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStep_Theta()
504 …PetscCall(TSComputeIHessianProductFunctionUP(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStep_Theta()
507 PetscCall(VecScale(VecsSensi2Temp[nadj], th->shift)); in TSAdjointStep_Theta()
524 if (th->endpoint) { /* two-stage Theta methods with th->Theta!=1, th->Theta==1 leads to BEuler */ in TSAdjointStep_Theta()
525 th->shift = 1. / ((th->Theta - 1.) * adjoint_time_step); in TSAdjointStep_Theta()
526 th->stage_time = adjoint_ptime; in TSAdjointStep_Theta()
527 PetscCall(TSComputeSNESJacobian(ts, th->X0, J, Jpre)); in TSAdjointStep_Theta()
530 if (quadts) PetscCall(TSComputeRHSJacobian(quadts, adjoint_ptime, th->X0, quadJ, NULL)); in TSAdjointStep_Theta()
540 PetscCall(VecScale(ts->vecs_sensi[nadj], 1. / th->shift)); in TSAdjointStep_Theta()
546 PetscCall(MatDenseGetColumn(th->MatFwdSensip0, 0, &xarr)); in TSAdjointStep_Theta()
549 …PetscCall(TSComputeIHessianProductFunctionUU(ts, adjoint_ptime, th->X0, VecsDeltaLam, ts->vec_sens… in TSAdjointStep_Theta()
551 PetscCall(MatDenseRestoreColumn(th->MatFwdSensip0, &xarr)); in TSAdjointStep_Theta()
553 …PetscCall(TSComputeIHessianProductFunctionUP(ts, adjoint_ptime, th->X0, VecsDeltaLam, ts->vec_dir,… in TSAdjointStep_Theta()
559 PetscCall(VecScale(ts->vecs_sensi2[nadj], 1. / th->shift)); in TSAdjointStep_Theta()
563 th->stage_time = ts->ptime; /* recover the old value */ in TSAdjointStep_Theta()
567 th->shift = 1.0 / (adjoint_time_step * th->Theta); in TSAdjointStep_Theta()
568 PetscCall(VecAXPBYPCZ(th->Xdot, -th->shift, th->shift, 0, th->X0, ts->vec_sol)); in TSAdjointStep_Theta()
569 …PetscCall(TSComputeIJacobianP(ts, th->stage_time, ts->vec_sol, th->Xdot, -1. / (th->Theta * adjoin… in TSAdjointStep_Theta()
570 if (quadts) PetscCall(TSComputeRHSJacobianP(quadts, th->stage_time, ts->vec_sol, quadJp)); in TSAdjointStep_Theta()
573 … PetscCall(VecAXPY(ts->vecs_sensip[nadj], -adjoint_time_step * th->Theta, VecsDeltaMu[nadj])); in TSAdjointStep_Theta()
577 … PetscCall(VecAXPY(ts->vecs_sensip[nadj], adjoint_time_step * th->Theta, ts->vec_drdp_col)); in TSAdjointStep_Theta()
587 …PetscCall(TSComputeIHessianProductFunctionPU(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStep_Theta()
592 …PetscCall(TSComputeIHessianProductFunctionPP(ts, th->stage_time, ts->vec_sol, VecsDeltaLam, ts->ve… in TSAdjointStep_Theta()
596 … PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * th->Theta, VecsDeltaMu2[nadj])); in TSAdjointStep_Theta()
597 …if (ts->vecs_fpu) PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * th->Theta, ts->ve… in TSAdjointStep_Theta()
598 …if (ts->vecs_fpp) PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * th->Theta, ts->ve… in TSAdjointStep_Theta()
603 PetscCall(VecZeroEntries(th->Xdot)); in TSAdjointStep_Theta()
604 …PetscCall(TSComputeIJacobianP(ts, adjoint_ptime, th->X0, th->Xdot, 1. / ((th->Theta - 1.0) * adjoi… in TSAdjointStep_Theta()
605 if (quadts) PetscCall(TSComputeRHSJacobianP(quadts, adjoint_ptime, th->X0, quadJp)); in TSAdjointStep_Theta()
608 …PetscCall(VecAXPY(ts->vecs_sensip[nadj], -adjoint_time_step * (1.0 - th->Theta), VecsDeltaMu[nadj]… in TSAdjointStep_Theta()
612 …PetscCall(VecAXPY(ts->vecs_sensip[nadj], adjoint_time_step * (1.0 - th->Theta), ts->vec_drdp_col)); in TSAdjointStep_Theta()
618 PetscCall(MatDenseGetColumn(th->MatFwdSensip0, 0, &xarr)); in TSAdjointStep_Theta()
621 …PetscCall(TSComputeIHessianProductFunctionPU(ts, adjoint_ptime, th->X0, VecsDeltaLam, ts->vec_sens… in TSAdjointStep_Theta()
623 PetscCall(MatDenseRestoreColumn(th->MatFwdSensip0, &xarr)); in TSAdjointStep_Theta()
625 …PetscCall(TSComputeIHessianProductFunctionPP(ts, adjoint_ptime, th->X0, VecsDeltaLam, ts->vec_dir,… in TSAdjointStep_Theta()
629 …PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * (1.0 - th->Theta), VecsDeltaMu2[nad… in TSAdjointStep_Theta()
630 …fpu) PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * (1.0 - th->Theta), ts->vecs_fp… in TSAdjointStep_Theta()
631 …fpp) PetscCall(VecAXPY(ts->vecs_sensi2p[nadj], -adjoint_time_step * (1.0 - th->Theta), ts->vecs_fp… in TSAdjointStep_Theta()
637 th->shift = 0.0; in TSAdjointStep_Theta()
638 PetscCall(TSComputeSNESJacobian(ts, th->X, J, Jpre)); /* get -f_y */ in TSAdjointStep_Theta()
640 if (quadts) PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, th->X, quadJ, NULL)); in TSAdjointStep_Theta()
653 th->shift = 1.0 / (adjoint_time_step * th->Theta); in TSAdjointStep_Theta()
654 PetscCall(VecAXPBYPCZ(th->Xdot, -th->shift, th->shift, 0, th->X0, th->X)); in TSAdjointStep_Theta()
655 …PetscCall(TSComputeIJacobianP(ts, th->stage_time, th->X, th->Xdot, th->shift, ts->Jacp, PETSC_FALS… in TSAdjointStep_Theta()
656 if (quadts) PetscCall(TSComputeRHSJacobianP(quadts, th->stage_time, th->X, quadJp)); in TSAdjointStep_Theta()
671 th->status = TS_STEP_COMPLETE; in TSAdjointStep_Theta()
677 TS_Theta *th = (TS_Theta *)ts->data; in TSInterpolate_Theta() local
681 PetscCall(VecCopy(ts->vec_sol, th->X)); in TSInterpolate_Theta()
682 if (th->endpoint) dt *= th->Theta; in TSInterpolate_Theta()
683 PetscCall(VecWAXPY(X, dt, th->Xdot, th->X)); in TSInterpolate_Theta()
689 TS_Theta *th = (TS_Theta *)ts->data; in TSEvaluateWLTE_Theta() local
691 Vec Y = th->vec_lte_work; /* Y = X + LTE */ in TSEvaluateWLTE_Theta()
695 if (!th->vec_sol_prev) { in TSEvaluateWLTE_Theta()
715 vecs[1] = th->X0; in TSEvaluateWLTE_Theta()
716 vecs[2] = th->vec_sol_prev; in TSEvaluateWLTE_Theta()
727 TS_Theta *th = (TS_Theta *)ts->data; in TSRollBack_Theta() local
731 if (quadts && ts->costintegralfwd) PetscCall(VecCopy(th->VecCostIntegral0, quadts->vec_sol)); in TSRollBack_Theta()
732 th->status = TS_STEP_INCOMPLETE; in TSRollBack_Theta()
733 if (ts->mat_sensip) PetscCall(MatCopy(th->MatFwdSensip0, ts->mat_sensip, SAME_NONZERO_PATTERN)); in TSRollBack_Theta()
734 …if (quadts && quadts->mat_sensip) PetscCall(MatCopy(th->MatIntegralSensip0, quadts->mat_sensip, SA… in TSRollBack_Theta()
740 TS_Theta *th = (TS_Theta *)ts->data; in TSForwardStep_Theta() local
742 Mat MatDeltaFwdSensip = th->MatDeltaFwdSensip; in TSForwardStep_Theta()
743 Vec VecDeltaFwdSensipCol = th->VecDeltaFwdSensipCol; in TSForwardStep_Theta()
751 previous_shift = th->shift; in TSForwardStep_Theta()
752 PetscCall(MatCopy(ts->mat_sensip, th->MatFwdSensip0, SAME_NONZERO_PATTERN)); in TSForwardStep_Theta()
754 …if (quadts && quadts->mat_sensip) PetscCall(MatCopy(quadts->mat_sensip, th->MatIntegralSensip0, SA… in TSForwardStep_Theta()
763 if (th->endpoint) { /* 2-stage method*/ in TSForwardStep_Theta()
764 th->shift = 1. / ((th->Theta - 1.) * th->time_step0); in TSForwardStep_Theta()
765 … PetscCall(TSComputeIJacobian(ts, th->ptime0, th->X0, th->Xdot, th->shift, J, Jpre, PETSC_FALSE)); in TSForwardStep_Theta()
767 PetscCall(MatScale(MatDeltaFwdSensip, (th->Theta - 1.) / th->Theta)); in TSForwardStep_Theta()
771 PetscCall(VecZeroEntries(th->Xdot)); in TSForwardStep_Theta()
772 …PetscCall(TSComputeIJacobianP(ts, th->ptime0, th->X0, th->Xdot, th->shift, ts->Jacp, PETSC_FALSE)); in TSForwardStep_Theta()
773 …PetscCall(MatAXPY(MatDeltaFwdSensip, (th->Theta - 1.) / th->Theta, ts->Jacp, SUBSET_NONZERO_PATTER… in TSForwardStep_Theta()
774 th->shift = previous_shift; in TSForwardStep_Theta()
775 PetscCall(VecAXPBYPCZ(th->Xdot, -th->shift, th->shift, 0, th->X0, ts->vec_sol)); in TSForwardStep_Theta()
776 …PetscCall(TSComputeIJacobianP(ts, th->stage_time, ts->vec_sol, th->Xdot, th->shift, ts->Jacp, PETS… in TSForwardStep_Theta()
780 th->shift = 0.0; in TSForwardStep_Theta()
781 …PetscCall(TSComputeIJacobian(ts, th->stage_time, th->X, th->Xdot, th->shift, J, Jpre, PETSC_FALSE)… in TSForwardStep_Theta()
787 th->shift = previous_shift; in TSForwardStep_Theta()
788 PetscCall(VecAXPBYPCZ(th->Xdot, -th->shift, th->shift, 0, th->X0, th->X)); in TSForwardStep_Theta()
789 …PetscCall(TSComputeIJacobianP(ts, th->stage_time, th->X, th->Xdot, th->shift, ts->Jacp, PETSC_FALS… in TSForwardStep_Theta()
795 th->shift = previous_shift; /* recover the previous shift used in TSStep_Theta() */ in TSForwardStep_Theta()
796 if (th->endpoint) { in TSForwardStep_Theta()
797 …PetscCall(TSComputeIJacobian(ts, th->stage_time, ts->vec_sol, th->Xdot, th->shift, J, Jpre, PETSC_… in TSForwardStep_Theta()
799 …PetscCall(TSComputeIJacobian(ts, th->stage_time, th->X, th->Xdot, th->shift, J, Jpre, PETSC_FALSE)… in TSForwardStep_Theta()
808 if (th->endpoint) { /* 2-stage method only */ in TSForwardStep_Theta()
810 PetscCall(TSComputeRHSJacobian(quadts, th->ptime0, th->X0, quadJ, NULL)); in TSForwardStep_Theta()
811 PetscCall(TSComputeRHSJacobianP(quadts, th->ptime0, th->X0, quadJp)); in TSForwardStep_Theta()
812 …PetscCall(MatTransposeMatMult(ts->mat_sensip, quadJ, MAT_REUSE_MATRIX, PETSC_DETERMINE, &th->MatIn… in TSForwardStep_Theta()
813 PetscCall(MatAXPY(th->MatIntegralSensipTemp, 1, quadJp, SAME_NONZERO_PATTERN)); in TSForwardStep_Theta()
814 …PetscCall(MatAXPY(quadts->mat_sensip, th->time_step0 * (1. - th->Theta), th->MatIntegralSensipTemp… in TSForwardStep_Theta()
819 for (ntlm = 0; ntlm < th->num_tlm; ntlm++) { in TSForwardStep_Theta()
823 if (th->endpoint) { in TSForwardStep_Theta()
846 if (!th->endpoint) { in TSForwardStep_Theta()
848 PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, th->X, quadJ, NULL)); in TSForwardStep_Theta()
849 PetscCall(TSComputeRHSJacobianP(quadts, th->stage_time, th->X, quadJp)); in TSForwardStep_Theta()
850 …PetscCall(MatTransposeMatMult(ts->mat_sensip, quadJ, MAT_REUSE_MATRIX, PETSC_DETERMINE, &th->MatIn… in TSForwardStep_Theta()
851 PetscCall(MatAXPY(th->MatIntegralSensipTemp, 1, quadJp, SAME_NONZERO_PATTERN)); in TSForwardStep_Theta()
852 …PetscCall(MatAXPY(quadts->mat_sensip, th->time_step0, th->MatIntegralSensipTemp, SAME_NONZERO_PATT… in TSForwardStep_Theta()
853 …PetscCall(MatAXPY(ts->mat_sensip, (1. - th->Theta) / th->Theta, MatDeltaFwdSensip, SAME_NONZERO_PA… in TSForwardStep_Theta()
855 PetscCall(TSComputeRHSJacobian(quadts, th->stage_time, ts->vec_sol, quadJ, NULL)); in TSForwardStep_Theta()
856 PetscCall(TSComputeRHSJacobianP(quadts, th->stage_time, ts->vec_sol, quadJp)); in TSForwardStep_Theta()
857 …PetscCall(MatTransposeMatMult(ts->mat_sensip, quadJ, MAT_REUSE_MATRIX, PETSC_DETERMINE, &th->MatIn… in TSForwardStep_Theta()
858 PetscCall(MatAXPY(th->MatIntegralSensipTemp, 1, quadJp, SAME_NONZERO_PATTERN)); in TSForwardStep_Theta()
859 …PetscCall(MatAXPY(quadts->mat_sensip, th->time_step0 * th->Theta, th->MatIntegralSensipTemp, SAME_… in TSForwardStep_Theta()
862 …if (!th->endpoint) PetscCall(MatAXPY(ts->mat_sensip, 1. / th->Theta, MatDeltaFwdSensip, SAME_NONZE… in TSForwardStep_Theta()
869 TS_Theta *th = (TS_Theta *)ts->data; in TSForwardGetStages_Theta() local
873 if (!th->endpoint && th->Theta != 1.0) *ns = 1; /* midpoint form */ in TSForwardGetStages_Theta()
877 if (!th->endpoint && th->Theta != 1.0) { in TSForwardGetStages_Theta()
878 th->MatFwdStages[0] = th->MatDeltaFwdSensip; in TSForwardGetStages_Theta()
880 th->MatFwdStages[0] = th->MatFwdSensip0; in TSForwardGetStages_Theta()
881 th->MatFwdStages[1] = ts->mat_sensip; /* stiffly accurate */ in TSForwardGetStages_Theta()
883 *stagesensip = th->MatFwdStages; in TSForwardGetStages_Theta()
890 TS_Theta *th = (TS_Theta *)ts->data; in TSReset_Theta() local
893 PetscCall(VecDestroy(&th->X)); in TSReset_Theta()
894 PetscCall(VecDestroy(&th->Xdot)); in TSReset_Theta()
895 PetscCall(VecDestroy(&th->X0)); in TSReset_Theta()
896 PetscCall(VecDestroy(&th->affine)); in TSReset_Theta()
898 PetscCall(VecDestroy(&th->vec_sol_prev)); in TSReset_Theta()
899 PetscCall(VecDestroy(&th->vec_lte_work)); in TSReset_Theta()
901 PetscCall(VecDestroy(&th->VecCostIntegral0)); in TSReset_Theta()
907 TS_Theta *th = (TS_Theta *)ts->data; in TSAdjointReset_Theta() local
910 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsDeltaLam)); in TSAdjointReset_Theta()
911 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsDeltaMu)); in TSAdjointReset_Theta()
912 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsDeltaLam2)); in TSAdjointReset_Theta()
913 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsDeltaMu2)); in TSAdjointReset_Theta()
914 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsSensiTemp)); in TSAdjointReset_Theta()
915 PetscCall(VecDestroyVecs(ts->numcost, &th->VecsSensi2Temp)); in TSAdjointReset_Theta()
945 TS_Theta *th = (TS_Theta *)ts->data; in SNESTSFormFunction_Theta() local
948 PetscReal shift = th->shift; in SNESTSFormFunction_Theta()
962 PetscCall(TSComputeIFunction(ts, th->stage_time, x, Xdot, y, PETSC_FALSE)); in SNESTSFormFunction_Theta()
970 TS_Theta *th = (TS_Theta *)ts->data; in SNESTSFormJacobian_Theta() local
973 PetscReal shift = th->shift; in SNESTSFormJacobian_Theta()
982 PetscCall(TSComputeIJacobian(ts, th->stage_time, x, Xdot, shift, A, B, PETSC_FALSE)); in SNESTSFormJacobian_Theta()
990 TS_Theta *th = (TS_Theta *)ts->data; in TSForwardSetUp_Theta() local
995 th->num_tlm = ts->num_parameters; in TSForwardSetUp_Theta()
996 PetscCall(MatDuplicate(ts->mat_sensip, MAT_DO_NOT_COPY_VALUES, &th->MatDeltaFwdSensip)); in TSForwardSetUp_Theta()
998 PetscCall(MatDuplicate(quadts->mat_sensip, MAT_DO_NOT_COPY_VALUES, &th->MatIntegralSensipTemp)); in TSForwardSetUp_Theta()
999 PetscCall(MatDuplicate(quadts->mat_sensip, MAT_DO_NOT_COPY_VALUES, &th->MatIntegralSensip0)); in TSForwardSetUp_Theta()
1002 PetscCall(MatDuplicate(ts->mat_sensip, MAT_DO_NOT_COPY_VALUES, &th->MatFwdSensip0)); in TSForwardSetUp_Theta()
1004 PetscCall(VecDuplicate(ts->vec_sol, &th->VecDeltaFwdSensipCol)); in TSForwardSetUp_Theta()
1010 TS_Theta *th = (TS_Theta *)ts->data; in TSForwardReset_Theta() local
1015 PetscCall(MatDestroy(&th->MatIntegralSensipTemp)); in TSForwardReset_Theta()
1016 PetscCall(MatDestroy(&th->MatIntegralSensip0)); in TSForwardReset_Theta()
1018 PetscCall(VecDestroy(&th->VecDeltaFwdSensipCol)); in TSForwardReset_Theta()
1019 PetscCall(MatDestroy(&th->MatDeltaFwdSensip)); in TSForwardReset_Theta()
1020 PetscCall(MatDestroy(&th->MatFwdSensip0)); in TSForwardReset_Theta()
1026 TS_Theta *th = (TS_Theta *)ts->data; in TSSetUp_Theta() local
1031 if (!th->VecCostIntegral0 && quadts && ts->costintegralfwd) { /* back up cost integral */ in TSSetUp_Theta()
1032 PetscCall(VecDuplicate(quadts->vec_sol, &th->VecCostIntegral0)); in TSSetUp_Theta()
1034 if (!th->X) PetscCall(VecDuplicate(ts->vec_sol, &th->X)); in TSSetUp_Theta()
1035 if (!th->Xdot) PetscCall(VecDuplicate(ts->vec_sol, &th->Xdot)); in TSSetUp_Theta()
1036 if (!th->X0) PetscCall(VecDuplicate(ts->vec_sol, &th->X0)); in TSSetUp_Theta()
1037 if (th->endpoint) PetscCall(VecDuplicate(ts->vec_sol, &th->affine)); in TSSetUp_Theta()
1039 th->order = (th->Theta == 0.5) ? 2 : 1; in TSSetUp_Theta()
1040 th->shift = 1 / (th->Theta * ts->time_step); in TSSetUp_Theta()
1050 if (!th->vec_sol_prev) PetscCall(VecDuplicate(ts->vec_sol, &th->vec_sol_prev)); in TSSetUp_Theta()
1051 if (!th->vec_lte_work) PetscCall(VecDuplicate(ts->vec_sol, &th->vec_lte_work)); in TSSetUp_Theta()
1055 ts->stifflyaccurate = (!th->endpoint && th->Theta != 1.0) ? PETSC_FALSE : PETSC_TRUE; in TSSetUp_Theta()
1061 TS_Theta *th = (TS_Theta *)ts->data; in TSAdjointSetUp_Theta() local
1064 PetscCall(VecDuplicateVecs(ts->vecs_sensi[0], ts->numcost, &th->VecsDeltaLam)); in TSAdjointSetUp_Theta()
1065 PetscCall(VecDuplicateVecs(ts->vecs_sensi[0], ts->numcost, &th->VecsSensiTemp)); in TSAdjointSetUp_Theta()
1066 …if (ts->vecs_sensip) PetscCall(VecDuplicateVecs(ts->vecs_sensip[0], ts->numcost, &th->VecsDeltaMu)… in TSAdjointSetUp_Theta()
1068 PetscCall(VecDuplicateVecs(ts->vecs_sensi[0], ts->numcost, &th->VecsDeltaLam2)); in TSAdjointSetUp_Theta()
1069 PetscCall(VecDuplicateVecs(ts->vecs_sensi2[0], ts->numcost, &th->VecsSensi2Temp)); in TSAdjointSetUp_Theta()
1075 PetscCall(VecDuplicateVecs(ts->vecs_sensi2p[0], ts->numcost, &th->VecsDeltaMu2)); in TSAdjointSetUp_Theta()
1085 TS_Theta *th = (TS_Theta *)ts->data; in TSSetFromOptions_Theta() local
1090 …s_theta_theta", "Location of stage (0<Theta<=1)", "TSThetaSetTheta", th->Theta, &th->Theta, NULL)); in TSSetFromOptions_Theta()
1091 …nstead of midpoint form of the Theta method", "TSThetaSetEndpoint", th->endpoint, &th->endpoint, N… in TSSetFromOptions_Theta()
1092 …previous solution (sometimes unstable)", "TSThetaSetExtrapolate", th->extrapolate, &th->extrapolat… in TSSetFromOptions_Theta()
1100 TS_Theta *th = (TS_Theta *)ts->data; in TSView_Theta() local
1106 PetscCall(PetscViewerASCIIPrintf(viewer, " Theta=%g\n", (double)th->Theta)); in TSView_Theta()
1107 … PetscCall(PetscViewerASCIIPrintf(viewer, " Extrapolation=%s\n", th->extrapolate ? "yes" : "no")); in TSView_Theta()
1114 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaGetTheta_Theta() local
1117 *theta = th->Theta; in TSThetaGetTheta_Theta()
1123 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaSetTheta_Theta() local
1127 th->Theta = theta; in TSThetaSetTheta_Theta()
1128 th->order = (th->Theta == 0.5) ? 2 : 1; in TSThetaSetTheta_Theta()
1134 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaGetEndpoint_Theta() local
1137 *endpoint = th->endpoint; in TSThetaGetEndpoint_Theta()
1143 TS_Theta *th = (TS_Theta *)ts->data; in TSThetaSetEndpoint_Theta() local
1146 th->endpoint = flg; in TSThetaSetEndpoint_Theta()
1154 TS_Theta *th = (TS_Theta *)ts->data; in TSComputeLinearStability_Theta() local
1157 f = (1.0 + (1.0 - th->Theta) * z) / (1.0 - th->Theta * z); in TSComputeLinearStability_Theta()
1166 TS_Theta *th = (TS_Theta *)ts->data; in TSGetStages_Theta() local
1170 if (!th->endpoint && th->Theta != 1.0) *ns = 1; /* midpoint form */ in TSGetStages_Theta()
1174 if (!th->endpoint && th->Theta != 1.0) { in TSGetStages_Theta()
1175 th->Stages[0] = th->X; in TSGetStages_Theta()
1177 th->Stages[0] = th->X0; in TSGetStages_Theta()
1178 th->Stages[1] = ts->vec_sol; /* stiffly accurate */ in TSGetStages_Theta()
1180 *Y = th->Stages; in TSGetStages_Theta()
1235 TS_Theta *th; in TSCreate_Theta() local
1269 PetscCall(PetscNew(&th)); in TSCreate_Theta()
1270 ts->data = (void *)th; in TSCreate_Theta()
1272 th->VecsDeltaLam = NULL; in TSCreate_Theta()
1273 th->VecsDeltaMu = NULL; in TSCreate_Theta()
1274 th->VecsSensiTemp = NULL; in TSCreate_Theta()
1275 th->VecsSensi2Temp = NULL; in TSCreate_Theta()
1277 th->extrapolate = PETSC_FALSE; in TSCreate_Theta()
1278 th->Theta = 0.5; in TSCreate_Theta()
1279 th->order = 2; in TSCreate_Theta()
1393 TS_Theta *th = (TS_Theta *)ts->data; in TSSetUp_BEuler() local
1396 …PetscCheck(th->Theta == 1.0, PetscObjectComm((PetscObject)ts), PETSC_ERR_OPT_OVERWRITE, "Can not c… in TSSetUp_BEuler()
1397 …PetscCheck(!th->endpoint, PetscObjectComm((PetscObject)ts), PETSC_ERR_OPT_OVERWRITE, "Can not chan… in TSSetUp_BEuler()
1431 TS_Theta *th = (TS_Theta *)ts->data; in TSSetUp_CN() local
1434 …PetscCheck(th->Theta == 0.5, PetscObjectComm((PetscObject)ts), PETSC_ERR_OPT_OVERWRITE, "Can not c… in TSSetUp_CN()
1435 …PetscCheck(th->endpoint, PetscObjectComm((PetscObject)ts), PETSC_ERR_OPT_OVERWRITE, "Can not chang… in TSSetUp_CN()