xcompact3d / x3d2

https://xcompact3d.github.io/x3d2
BSD 3-Clause "New" or "Revised" License
3 stars 4 forks source link

Add Adams-Bashforth time-integration schemes #88

Closed CFD-Xing closed 3 months ago

CFD-Xing commented 3 months ago

This PR implements Adams-Bashforth integration schemes with provision for the start-up phase by successively increasing the time-integration order.

The Dahlquist equation is used to check Adams-Bashforth integration schemes from order 1 to 4. New tests have been added in the following file tests/omp/test_omp_adamsbashforth.f90

CFD-Xing commented 3 months ago

Looking good. It would be good to have tests for each time advancement scheme (maybe using a simple equation with an analytical solution, to check that they indeed converge to the expected order). RK schemes (at least 3rd order and possibly 4th order) would be good to have, possibly in another PR.

Implementing RK schemes is going to be more challenging and will require more significant changes in order to work with a segregated solver. The velocity update, Poisson solver, and velocity correction need to be solved in stages, this will require changes to both solver.f90 and time_integration.f90. Basically, the whole time-integration process need to be wrapped into a function/subroutine.

semi-h commented 3 months ago

@Nanoseb also mentioned earlier but it looks like there is a pattern across adams_bashforth 1 to 4 and this can be handled with a single subroutine if we pass the order parameter. The coefficients belong to different orders can be stored in an array (or array of pointers) and then a loop over them can sort the vecadd lines, and then in the part the pointers are rolled the if check is based on the order variable anyways so it looks quite straightforward to implement.

Also, Fortran indexing is 1 based so starting istep from 1 could be a better idea, any thoughts on this?

pbartholomew08 commented 3 months ago

Also, Fortran indexing is 1 based so starting istep from 1 could be a better idea, any thoughts on this?

I could go either way - we aren't indexing into an array with istep (e.g. a record of all timesteps) and time starts at t=0, so in my opinion istep=0 makes sense (also gives you t = istep * dt). However in Xcompact3d istep starts at 1 so from familiarity perspective it might be better to retain this.

slaizet commented 3 months ago

I would prefer to start at 1 please :-)

Nanoseb commented 3 months ago

see failed test. Note that you can run them locally with make test.

CFD-Xing commented 3 months ago

enstrophy

This is the enstrophy figure obtained with this PR