Skip to content

Commit e740fc5

Browse files
committed
Land the last time step exactly on the final time
The transient loop used `for(T=0; T<=FinalTime;)` with `T+=Dt` at the end of a successful step, so it always ran one step too many and finished at up to FinalTime+Dt instead of FinalTime (e.g. end-time=0.01 with dt=0.003 finished at 0.012, a 20% overshoot; with adaptive dt the overshoot can be a full dtmax). Use a strict bound with a small relative tolerance and clamp the increment of the final step so it lands exactly on the final time. Verified: a fixed-dt run now ends at 1.00000e-02 (was 1.20000e-02); adaptive diffusion and J2-plasticity runs land on their end-time; the clamped result is partition-consistent (np=1 vs np=4 agree to 6e-17).
1 parent 73b5c23 commit e740fc5

1 file changed

Lines changed: 8 additions & 1 deletion

File tree

src/TimeStepping/TimeSteppingSolve.cpp

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -80,7 +80,14 @@ bool TimeStepping::solve(FECell &t_FECell,
8080
MessagePrinter::printStars();
8181

8282
IsLastStepFailed=true;
83-
for(t_FECtrlInfo.T=0.0;t_FECtrlInfo.T<=m_Data.m_FinalTime;){
83+
// tolerance (relative to the final time) for deciding we have reached the final time
84+
const double finaltimetol=1.0e-10*m_Data.m_FinalTime;
85+
for(t_FECtrlInfo.T=0.0;t_FECtrlInfo.T<m_Data.m_FinalTime-finaltimetol;){
86+
// clamp the increment so the last step lands exactly on the final time instead of
87+
// overshooting it by up to one dt
88+
if(t_FECtrlInfo.T+t_FECtrlInfo.Dt>m_Data.m_FinalTime){
89+
t_FECtrlInfo.Dt=m_Data.m_FinalTime-t_FECtrlInfo.T;
90+
}
8491
snprintf(buff,68,"Time=%13.5e, step=%8d, dt=%13.5e",t_FECtrlInfo.T+t_FECtrlInfo.Dt,t_FECtrlInfo.CurrentStep+1,t_FECtrlInfo.Dt);
8592
str=buff;
8693
MessagePrinter::printNormalTxt(str);

0 commit comments

Comments
 (0)