Two Adaptive Timestepping questions/suggestions
-
- Posts: 196
- Joined: 29 Sep 2011, 12:25
- Antispam: Yes
Two Adaptive Timestepping questions/suggestions
Hello Everbody,
Hello dear Elmer Team,
i hope everybody had a good summer and that you guys missed me as much as i missed you. i had a bit of an Elmer braeak the last months but will start hacking around slowly now.
I just read through the ElmerSolver file and the adaptive timestepping therein. I am a big fan of thew "two small steps, one big , compare" scheme, that i used myself already.
Ihave two questions to the Elmer implementation now:
1)
When the step is accepted wouldn't it be better to take the solution from the second small step? because this comes from a finer time discretisation?
I think that would give some cheap performance improvement in transient simulations.
2)
The second thing is the way the error is calculated. it is calculated by the difference of the norms. Wouldn't the norm of the differences be a better measure? But maybe it's hard to calculate the norm of the difference in the Elmersolver module?
As always please don't blame me if i got everything wrong.
i am looking forward to your reply and wish a nice weekend.
best regards
Franz
Hello dear Elmer Team,
i hope everybody had a good summer and that you guys missed me as much as i missed you. i had a bit of an Elmer braeak the last months but will start hacking around slowly now.
I just read through the ElmerSolver file and the adaptive timestepping therein. I am a big fan of thew "two small steps, one big , compare" scheme, that i used myself already.
Ihave two questions to the Elmer implementation now:
1)
When the step is accepted wouldn't it be better to take the solution from the second small step? because this comes from a finer time discretisation?
I think that would give some cheap performance improvement in transient simulations.
2)
The second thing is the way the error is calculated. it is calculated by the difference of the norms. Wouldn't the norm of the differences be a better measure? But maybe it's hard to calculate the norm of the difference in the Elmersolver module?
As always please don't blame me if i got everything wrong.
i am looking forward to your reply and wish a nice weekend.
best regards
Franz
-
- Posts: 196
- Joined: 29 Sep 2011, 12:25
- Antispam: Yes
Re: Two Adaptive Timestepping questions/suggestions
OOops, SOrry.
Forget about my first point. Elmer does it the good way, of course! Sorry for being a wise guy.
Still the second would be interesting.
And as another suggestion
3)
what about introducing another limit. so that if the error is between the two limits the step is accepted but the timestep is not changed.
if the error is smaller than both limits the step is accepted and the timestep is increased.
and if the error is bigger than both limits than step is not accepted and timestep is decreased.
That way if the limits are chosen carefully often the simulation "finds" a good timestep by itself.
so
sorry for the mistake again
and have a nice one
Franz
Forget about my first point. Elmer does it the good way, of course! Sorry for being a wise guy.
Still the second would be interesting.
And as another suggestion
3)
what about introducing another limit. so that if the error is between the two limits the step is accepted but the timestep is not changed.
if the error is smaller than both limits the step is accepted and the timestep is increased.
and if the error is bigger than both limits than step is not accepted and timestep is decreased.
That way if the limits are chosen carefully often the simulation "finds" a good timestep by itself.
so
sorry for the mistake again
and have a nice one
Franz
Re: Two Adaptive Timestepping questions/suggestions
Hi,
the problem with the adaptive timestepping as it is now in Elmer is that it often makes the simulation time explode. As in some kinds of simulations, small timesteps are needed in order to achieve convergence at all, a suggestion from my side is to take convergence as additional criterion: Decrease timestep if solution did not converge.
I have implemented that one in a custom solver and use it frequently.
Matthias
the problem with the adaptive timestepping as it is now in Elmer is that it often makes the simulation time explode. As in some kinds of simulations, small timesteps are needed in order to achieve convergence at all, a suggestion from my side is to take convergence as additional criterion: Decrease timestep if solution did not converge.
I have implemented that one in a custom solver and use it frequently.
Matthias
Re: Two Adaptive Timestepping questions/suggestions
@mzenker: Could you share your code? The current adaptive timestepping continues to bisect the time until it reaches the limit in my cases (MEMS). I was just getting ready to try to write an adaptive timestepper dependent on convergence and seeing your code would really be great.
Re: Two Adaptive Timestepping questions/suggestions
Hi,
my solver does actually some other things as well, but I can share the relevant code pieces.
The idea is to check convergence after each timestep. If the calculation did not converge, the time and result is reset to previous, and the timestep is divided by two. If Resultoutput solver is used, a variable is set in order not to save the nonconverged result. As I don't use ElmerPost any more, I didn't take care of .ep output.
All this would ideally be done in main Elmersolver, but as I was not able to compile it at the time, I implemented the functionality in a custom solver.
I have made an extract of my solver which should only do the adaptive timestepping. I did not test it, so there may be some bugs and/or syntax errors in it. Here it is:
In ResultoutputSolver, I have added a piece of code which gets the SaveResults variable and returns without saving if it is zero:
Now in the sif file, the following is needed:
I think that's all, but I have not tested it. My actual implementation which works is almost identical.
I hope it will be useful.
Question to the Elmer team: What do you think about implementing this functionality as an extension to the existing adaptive timestepping in Elmersolver?
Matthias
EDIT (7 Nov 2012): Removed some confusing remainders from full solver, added initialization for "CurrentTimeStepSize"
my solver does actually some other things as well, but I can share the relevant code pieces.
The idea is to check convergence after each timestep. If the calculation did not converge, the time and result is reset to previous, and the timestep is divided by two. If Resultoutput solver is used, a variable is set in order not to save the nonconverged result. As I don't use ElmerPost any more, I didn't take care of .ep output.
All this would ideally be done in main Elmersolver, but as I was not able to compile it at the time, I implemented the functionality in a custom solver.
I have made an extract of my solver which should only do the adaptive timestepping. I did not test it, so there may be some bugs and/or syntax errors in it. Here it is:
Code: Select all
SUBROUTINE AdaptiveSolver(Model, Solver, dt, Transient)
USE DefUtils
IMPLICIT NONE
TYPE(Solver_t) :: Solver
TYPE(Model_t) :: Model
REAL(KIND=dp) :: dt
LOGICAL :: Transient
!----------------------------------------------------------------
! Local variables
!----------------------------------------------------------------
TYPE(ValueList_t), POINTER :: SolverParams
TYPE(Mesh_t), POINTER :: Mesh
TYPE(Variable_t), POINTER :: Var
REAL(KIND=dp), POINTER :: WrkPntr(:)
INTEGER :: SaveResults = 1
REAL(KIND=dp) :: MaxTimeStepSize = 0.0_dp, MinTimeStepSize, CurrentTimeStepSize=0.0_dp, MinSaveInterval
REAL(KIND=dp) :: Time, Last_Time = 0.0_dp, Last_Save_Time = 0.0_dp, Duration
REAL(KIND=dp) :: SaveResults_real
LOGICAL :: Found, Visited = .FALSE.
LOGICAL :: AdaptiveTimeStepping, Converged, StopSimulation = .FALSE.
INTEGER :: istat
TYPE(Solver_t), POINTER :: iSolver
REAL(KIND=dp), ALLOCATABLE :: xx(:,:), prevxx(:,:,:)
INTEGER :: i, j, k, n
INTEGER :: AdaptiveKeepSmallest, KeepTimeStepSize = 0
INTEGER :: SteadyConverged
SAVE Last_Time
SAVE CurrentTimeStepSize
SAVE KeepTimeStepSize
SAVE Last_Save_Time
SAVE Visited
SAVE xx, prevxx
SAVE StopSimulation
CALL Info( 'AdaptiveSolver', '*************************************************' )
!Read in solver parameters
!-------------------------
SolverParams => GetSolverParams()
IF (.NOT. ASSOCIATED(SolverParams)) CALL FATAL('AdaptiveSolver','No Solver section found')
IF ( .NOT. Visited ) THEN
! allocate space for previous result
n = Model % NumberOfSolvers
j = 0
k = 0
DO i=1,n
iSolver => Model % Solvers(i)
IF ( ASSOCIATED( iSolver % Variable % Values ) ) THEN
IF ( ASSOCIATED( iSolver % Variable % PrevValues ) ) THEN
j = MAX( j, SIZE( iSolver % Variable % PrevValues,2 ) )
END IF
k = MAX( k, SIZE( iSolver % Variable % Values ) )
END IF
END DO
ALLOCATE( xx(n,k), prevxx( n,k,j ) )
! initialize variables, allocate memory
Mesh => Model % Meshes
DO WHILE( ASSOCIATED(Mesh) )
IF ( Mesh % OutputActive ) EXIT
Mesh => Mesh % Next
END DO
CALL SetCurrentMesh( Model, Mesh )
NULLIFY(WrkPntr)
ALLOCATE(WrkPntr(1),STAT=istat)
IF( istat /= 0 ) CALL Fatal('AdaptiveSolver','Memory allocation error')
CALL VariableAdd( Model % Variables, Mesh, Solver, 'AdaptiveTimeStepSize', 1, WrkPntr )
NULLIFY(WrkPntr)
ALLOCATE(WrkPntr(1),STAT=istat)
IF( istat /= 0 ) CALL Fatal('AdaptiveSolver','Memory allocation error 3')
CALL VariableAdd( Model % Variables, Mesh, Solver, 'SaveResults', 1, WrkPntr )
END IF
Time = GetTime()
IF ( StopSimulation ) THEN
! Set SaveResults variable to zero
Var => VariableGet( Model % Variables, 'SaveResults' )
IF(.NOT. ASSOCIATED(Var)) THEN
CALL FATAL( 'AdaptiveSolver', 'No memory allocated for SaveResults')
END IF
! It seems impossible to save an Integer, so we convert it to Real
Var % Values(1) = 0.0_dp
! End of simulation was reached at previous timestep - now send stop signal
WRITE(Message, *) 'End of Simulation! Sending stop signal to solver...'
CALL ListAddConstReal(Model % Simulation,'Exit Condition',1.0_dp)
ELSE
Converged = .TRUE.
SaveResults = 1
DO i=1,Model % NumberOfSolvers
iSolver => Model % Solvers(i)
SteadyConverged = iSolver % Variable % SteadyConverged
IF( ( Time > 0 ) .AND. ( SteadyConverged == 0 ) ) THEN
Converged = .FALSE.
CALL Info('AdaptiveSolver', 'No convergence for this timestep')
EXIT
END IF
END DO
AdaptiveTimeStepping = ListGetLogical(SolverParams,'Convergence Adaptive Timestepping',Found)
IF ( AdaptiveTimeStepping ) THEN
MinTimestepSize = ListGetConstReal( SolverParams, 'Convergence Adaptive Min Timestep', Found )
IF ( .NOT. Found ) MinTimeStepSize = 1e-15
AdaptiveKeepSmallest = GetInteger(SolverParams, 'Convergence Adaptive Keep Smallest', Found )
IF ( .NOT. Found ) AdaptiveKeepSmallest = 2
IF( Converged .OR. CurrentTimeStepSize - MinTimeStepSize < 1e-10 ) THEN
! save current result
DO i=1,Model % NumberOFSolvers
iSolver => CurrentModel % Solvers(i)
IF ( ASSOCIATED( iSolver % Variable % Values ) ) THEN
n = SIZE( iSolver % Variable % Values )
xx(i,1:n) = iSolver % Variable % Values
IF ( ASSOCIATED( iSolver % Variable % PrevValues ) ) THEN
DO j=1,SIZE( iSolver % Variable % PrevValues,2 )
prevxx(i,1:n,j) = iSolver % Variable % PrevValues(:,j)
END DO
END IF
END IF
END DO
ELSE
! - reset to result for last convergence
DO i=1,Model % NumberOFSolvers
iSolver => CurrentModel % Solvers(i)
IF ( ASSOCIATED( iSolver % Variable % Values ) ) THEN
n = SIZE(iSolver % Variable % Values)
iSolver % Variable % Values = xx(i,1:n)
IF ( ASSOCIATED( iSolver % Variable % PrevValues ) ) THEN
DO j=1,SIZE( iSolver % Variable % PrevValues,2 )
iSolver % Variable % PrevValues(:,j) = prevxx(i,1:n,j)
END DO
END IF
END IF
END DO
! - reset time to previous value
Var => VariableGet( Model % Solver % Mesh % Variables, 'Time' )
IF ( ASSOCIATED( Var ) ) Var % Values(1) = Time - CurrentTimeStepSize
Time = GetTime()
Write(Message, *) 'Adaptive Timestepping: Time reset to ', Time
CALL Info('AdaptiveSolver', Message)
SaveResults = 0
END IF
END IF ! AdaptiveTimeStepping
! Get duration
Duration = ListGetConstReal( SolverParams, 'Duration',Found)
IF ( .NOT. Found ) Duration = 0.0
! Get time minimum saving interval
MinSaveInterval = ListGetConstReal(SolverParams, 'Min Data Saving Interval',Found)
IF ( .NOT. Found ) MinSaveInterval = 0.0_dp
IF( SaveResults /= 0 ) THEN
IF( ( Time > 1e-10 ) .AND. ( Time - Last_Save_Time < MinSaveInterval ) ) THEN
SaveResults = 0
ELSE
Last_Save_Time = Time
END IF
END IF
! Check if duration is exceeded since Time was saved the last time
! Duration is multiplied by a factor to avoid problems due to rounding errors
IF ( Time - Last_Time >= Duration - 1e-10 ) THEN
! End of simulation after this timestep! Save one last time...
SaveResults = 1
StopSimulation = .TRUE.
END IF
! Get time step size
MaxTimeStepSize = ListGetConstReal(SolverParams, 'Adaptive TimeStep Size',Found)
IF ( .NOT. Found ) MaxTimeStepSize = 0.0
! Determine timestep size
IF ( .NOT. Visited ) THEN
CurrentTimeStepSize = MaxTimeStepSize
ELSE IF ( AdaptiveTimeStepping ) THEN
IF ( Converged ) THEN
! check keepsmallest, increase timestep size if desired
IF ( KeepTimeStepSize > 0 ) THEN
KeepTimeStepSize = KeepTimeStepSize - 1
ELSE
! try to double the timestep until the timestep set by the user is reached
CurrentTimeStepSize = CurrentTimeStepSize * 2
IF ( CurrentTimeStepSize > MaxTimeStepSize ) THEN
CurrentTimeStepSize = MaxTimeStepSize
ELSE
WRITE( Message, * ), 'Adaptive Timestepping: Time step set to ', CurrentTimeStepSize
CALL Info('AdaptiveSolver', Message)
END IF
END IF
ELSE
! set timestep size to half if MinTimestepSize is not reached
IF( CurrentTimeStepSize - MinTimeStepSize > 1e-10 ) THEN
CurrentTimeStepSize = CurrentTimeStepSize / 2
IF ( CurrentTimeStepSize < MinTimeStepSize ) CurrentTimeStepSize = MinTimeStepSize
WRITE( Message, * ), 'Adaptive Timestepping: Time step set to ', CurrentTimeStepSize
CALL Info('AdaptiveSolver', Message)
ELSE
CALL Warn('AdaptiveSolver','Cannot set smaller time step!')
CALL Warn('AdaptiveSolver','Minimum time step might be too big to reach convergence.')
CALL Warn('AdaptiveSolver','Continuing anyway...')
END IF
! reset KeepTimeStepSize
KeepTimeStepSize = AdaptiveKeepSmallest
END IF
END IF
! save TimeStepSize
Var => VariableGet( Model % Variables, 'AdaptiveTimeStepSize' )
IF(.NOT. ASSOCIATED(Var)) CALL FATAL( 'AdaptiveSolver', 'No memory allocated for AdaptiveTimeStepSize')
Var % Values(1) = CurrentTimeStepSize
! Save SaveResults variable
Var => VariableGet( Model % Variables, 'SaveResults' )
IF(.NOT. ASSOCIATED(Var)) CALL FATAL( 'AdaptiveSolver', 'No memory allocated for SaveResults')
! It seems impossible to save an Integer, so we convert it to Real
SaveResults_real = REAL(SaveResults)
Var % Values(1) = SaveResults_real
IF ( .NOT. StopSimulation ) THEN
WRITE(Message, *) ' SaveResults', SaveResults
ELSE
WRITE(Message, *) 'This was the last timestep! Preparing end of Simulation...'
END IF
END IF
CALL Info( 'AdaptiveSolver', Message, level=4)
CALL Info( 'AdaptiveSolver', '*************************************************' )
IF ( .NOT. Visited ) Visited = .TRUE.
END SUBROUTINE AdaptiveSolver
Code: Select all
! Check SaveResults variable from Adaptive solver, don't save if it is equal to 0
IF( ListGetLogical(Params,'Check SaveResults Variable',Found) ) THEN
Var => VariableGet( Model % Variables, 'SaveResults' )
IF(.NOT. ASSOCIATED(Var)) THEN
CALL Warn( 'ResultOutputSolver', 'SaveResults not found' )
SaveResults = 1.0_dp
ELSE
SaveResults = Var % Values(1)
END IF
IF( SaveResults < 0.1 ) THEN
CALL Info( 'ResultOutputSolver', 'Nothing to save...')
CALL Info( 'ResultOutputSolver', '-------------------------------------')
RETURN
END IF
END IF
Code: Select all
Simulation
...
Timestep intervals = 10000
Timestep Size = Equals AdaptiveTimeStepSize
...
End
...
Solver 1
Equation = Adaptive Solver
Convergence Adaptive Timestepping = Logical True
Convergence Adaptive Min Timestep = Real 0.1
Convergence Adaptive Keep Smallest = Integer 5
Duration = Real 1200
Variable = -dofs 1 -global AdaptiveTimeStepSize
Procedure = "AdaptiveSolver" "AdaptiveSolver"
Adaptive TimeStep Size = 5
Min Data Saving Interval = Real 1
Exec Solver = After Timestep
End
Solver 2
Equation = Result Output
Procedure = "ResultOutputSolve" "ResultOutputSolver"
...
Check SaveResults Variable = Logical True
Exec Solver = After Timestep
End
I hope it will be useful.
Question to the Elmer team: What do you think about implementing this functionality as an extension to the existing adaptive timestepping in Elmersolver?
Matthias
EDIT (7 Nov 2012): Removed some confusing remainders from full solver, added initialization for "CurrentTimeStepSize"
Last edited by mzenker on 08 Nov 2012, 17:58, edited 2 times in total.
Re: Two Adaptive Timestepping questions/suggestions
Wow, thanks! It would be excellent if this functionality were integrated into ElmerSolver... in the mean time, I did a compile check and got the following:
Thank you very much for sharing. I am starting to look at your code to try to understand it. Are IncState, Loop, and State variables part of the adaptive timestepper or part of your other code?
Code: Select all
IF ( IncState .OR. ( .NOT. Visited ) ) THEN
1
Error: Symbol 'incstate' at (1) has no IMPLICIT type
mzenkersolver.f90:229.53:
WRITE(Message, *) 'Next Step: Loop ', Loop, ' State ',State, ' SaveR
1
Error: Symbol 'loop' at (1) has no IMPLICIT type
mzenkersolver.f90:229.70:
WRITE(Message, *) 'Next Step: Loop ', Loop, ' State ',State, ' SaveR
1
Error: Symbol 'state' at (1) has no IMPLICIT type
Re: Two Adaptive Timestepping questions/suggestions
Hi,
IncState, State and Loop are part of the other code. You can replace
by
and remove Loop and State from the WRITE statements (and also from the comments).
Matthias
IncState, State and Loop are part of the other code. You can replace
Code: Select all
IF ( IncState .OR. ( .NOT. Visited ) ) THEN
Code: Select all
IF ( .NOT. Visited ) THEN
Matthias
-
- Posts: 196
- Joined: 29 Sep 2011, 12:25
- Antispam: Yes
Re: Two Adaptive Timestepping questions/suggestions
Hi there,
i would like to stress my second question once again, hoping to start a discussion about that.
The way the error is calculated.
Right now Elmer is taking the norms of the solutions and taking the difference between them as the error measure.
Wouldn't it be better to use the norm of the difference?
best regards Franz
i would like to stress my second question once again, hoping to start a discussion about that.
The way the error is calculated.
Right now Elmer is taking the norms of the solutions and taking the difference between them as the error measure.
Wouldn't it be better to use the norm of the difference?
best regards Franz
-
- Posts: 196
- Joined: 29 Sep 2011, 12:25
- Antispam: Yes
Re: Two Adaptive Timestepping questions/suggestions
Hi, this one is for mzenker,
i had a look at your solver from abovew. looks really nice.
i think i will use it too, so thank you for sharing it.
just a short remark:
in line 114 you use CurrentTimeStepSize uninitialized in the first run, i think.
could that be true?
mybe your full solver hasn't gotr that problem.
best regards
Franz
i had a look at your solver from abovew. looks really nice.
i think i will use it too, so thank you for sharing it.
just a short remark:
in line 114 you use CurrentTimeStepSize uninitialized in the first run, i think.
could that be true?
mybe your full solver hasn't gotr that problem.
best regards
Franz
Re: Two Adaptive Timestepping questions/suggestions
Hi Franz,
you are right, I have corrected that in the code. The missing initialization is also in my full solver and has not been a problem so far, maybe Fortran automatically initializes to zero. But it is safer not to rely to that, of course...
Thank you for pointing that out!
Matthias
you are right, I have corrected that in the code. The missing initialization is also in my full solver and has not been a problem so far, maybe Fortran automatically initializes to zero. But it is safer not to rely to that, of course...
Thank you for pointing that out!
Matthias