jothurgood/simplecfd

Simple (and not-so-simple) CFD solvers written in Fortran with Python plotting routines

simplecfd's Issues

CCAPS (all): Refactor so that each problem setup/diagnostic can be self contained in own module

Architecture for pre-defined test problems has become a bit unwieldy.

Would like to within the main code essentially have user control / initial conditions, but a one line override that points to a specific problem module that has its own control, initial condition and diagnostic that takes precedent over the main module.

CCAPS - use time centred rho in step 4's grav forcing

In step 4 should probabally calc a time centred rho (from 0.5 * (rho^n + rho^n+1 ) to use in calc of U*

Currently uses rho^n+1 (as advect_dens is called prior, which updates the rho array to the new time level).

Should test and then (if good) implement in CCAPS, CCAPS_vardens and CCAPS_atm_const

CCAPS - test "slip" bc and roll out to all vars

currently not properly tested or set up for all faces and flavours.

advection + burgers - bug if xmin neq 0?

check this and fix if need be

make relaxation step in CCAPS more modular

make relaxation step in CCAPS more modular instead of two calls to very similar subroutine. Then the simple projection tests can be built from the same makefile in a cleaner manner, and a future multigrid solver could just be parachuted in-line into the advance_module

enforce di=1 -> x and d=2 >y throughout the various functions consistently

some use di=0 for x components, d=1 for y. Others use 1 and 2. Fortran would typically use 1 and 2. Switch all to 1 and 2

MG - extend to full V cycle

add const gravitational forcing to CCAPS so that RTI test can be performed

want to test vardens on a Rayleigh Taylor instability

enforce a column width of 72

mainly because I like to have a terminal vertically split on a small laptop with fairly large font, so this reads well and doesn't visually break lines

CCAPS_vardens - incorporate new bc and problems

• update for the new boundary handling in CCAPS
• copy the circular drop problem in and check identical

MG - add time-independent test for the const coeff helmholtz

Current test is time-stepping a diffusion equation and doesnt give 2nd order convg.

Devise an independent test to check the order of MG itself when running in this mode is still 2nd.

CCAPS_atm_const - possible correction bug

Not seen any bad behaviour but compiler debug warning

[ 50%] Building Fortran object CMakeFiles/CCAPS.dir/src/advance_module.f90.o
/Users/tflv6/Data/SimpleCFD/incompressible/CCAPS_atm_const/src/advance_module.f90:313:0:

         betay_r = 0.5_num * (betacc + p0(iy+1)**(1.0_num / gamma))

Warning: 'betacc' may be used uninitialized in this function [-Wmaybe-uninitialized]

Inspecting code, seems will in fact occur when evaluating the iy=0 slice for ix>0.

This might be causing problems, although you've not spotted anything so far.

MG - test3 - is it okay?

test 3 is test 2 but passed to var eta solver for eta = 1 everywhere

eta bcs are set to "fixed" but fixed eta for inhomo not implemented. I.e., I think it might be seeing eta=0 on the boundaries.

Why does it still work and not cause problems?

remove stray debug code from ccaps advance_module

Some weird debugging condition on the output has crept into the main source, clean it up

!DEBUG
if (step >= 1211) then
    print *, '*** max divu before cleaning',maxval(abs(divu)*dt)
endif

Write an exact Riemann solver for 1D Euler equations in an ideal gas

Write an exact R solver for the 1D Euler equations in an ideal gas

To provide exactly converged solutions against which to compare the various numerical solutions.

(later ... possibly make a split 2D extension, or even MHD if feeling brave / have time)

Add all pedagogical CFD codes you have written

Add one at a time, from simple advection up, refactoring as go so all are consistent with one another

CCAPS - double check the gravitational dt

MG - generalise Poission to general elliptic

MG - make basic properly modular and incorporate into CCAPS

CCAPS - avoid overspecification of boundaries

The need to have vface bcs is just related to applying minmod at the very edge, otherwise the two ghosts for the sln vars should be enough to directly give you everything you need for the edge vars.

Add handling of the extreme case to minmod (just revert to a simple gradient there), and then change some of the loop sizes to capture the limits appropriately.

If you run with appropriate compiler flags, should check if OOB for you.

The simplification of boundary conditions should remove a number of possible sources of error/bugs and make the code easier to use.

Maths in wiki has stopped rendering properly

Problem with website that is meant to generate the html?

CCAPS (all flavours) - RTI bug (and avoid hard-coding of fixed eta for MG boundary)

Hard-coding fixed eta values for the MG boundary will be a hassle in the long run.

Also, current hardcoding corresponding to rho = 1 and rho = 7 for RTI isn't really correct? it should be those ^-1 . Why does the RTI mostly work? Why did you find you start to get nonconvergence issues when you correct it?

CCAPS - make all loops column-major

mixed order of access

Double check that you havent missed any 1/rho type factors any where in the algorithm

Refactor all -> remove pyplot-fortran in favour of directly calling own scripts

When first wrote these was using

https://github.com/jacobwilliams/pyplot-fortran

Useful, but gives it a bit of a modular dependency and that FoBis nonsense...

In interests of portability of the whole repository, better to just write own python scripts to read *dat and call from within Fortran so that images are automatically produced during run time (if so desired)

CCAPS : noslip + grav issue

in grav case, you still get some bad divergence left across the "opposite" no slip boundary to the direction of gravity (e.g. grav_y < 0 => y_max (no_slip) having this issue).

CCAPS - possible missing rho bc application in step_2

In advance_module

336 if (use_vardens) call rho_bcs ! needed for any OOB in relax and correction

but rho_bc was updated to need arguments about whether cc, x_face, y_face, etc.

Is this call just not updating the boundary?

Could be a big source of error!!!

CCAPS (all flavours) - check all factors of rho in code and wiki

CCAPS (all flavours) - add reporting of MG boundary conditions to setup_report

To keep track of options used in a run

CCAPS - remove (redundant) warning messages that MG is hardcoded to periodic

no longer true

CCAPS - No Slip + Grav issue

Description

• Set the RTI problem up with no density interface (e.g. rho = 1 throughout)
• Will see divergence + edge-mode v velocity emanate from top boundary (boundary opposite to grav drn)
• from printouts, correction to v is not correctly implemented near top boundary (at iy = ny-1 and iy = ny, but with the cells at iy=ny being a much larger error)

Cause

• the smaller errors at iy = ny-1 are a consequence of the BC choice of extrapolation of pressure through the boundary. Different BC choices / higher order interpolation / not using centred differences for evaluation of pressure gradient near the edges might improve this (I suspect not using CD but "inward difference" near the edge is the way to go so gradp is never determined based on boundary?)

• the larger error at iy = ny is mainly due to rho(ix,iy=ny) not remaining =1 at the boundary (although the choice of pressure bc also effects it somewhat)

• this is an unexpected result from advect_dens

• Tracing it back, it emanates from the evaluation of the right (full/transversal) rho state)

705         rhoyr(ix,iy) = rhohyr(ix,iy) &
706           & - 0.5_num * dt * rho(ix,iy+1) * (macv(ix,iy+1) - macv(ix,iy)) / dy &
707           & - 0.5_num * dt / dx * ( rhohx(ix,iy+1)*macu(ix,iy+1) - &
708                                     & rhohx(ix-1,iy+1)*macu(ix-1,iy+1) )  

• Specifically the term at 706.

• This is because you "expected" macv(ix,iy+1) - macv(ix,iy) to evaluate as 0 but no slip means that macv(iy+1) = - macv(iy), so you get a non zero term

• So the boundary as is doesn't really make sense with the gravity? Could either implement a zero-gradient / inoutflow boundary at the top, or implement (after this step, gradp should correctly compensate the gravity and noslip on internal flows makes more sense?)

CCAPS IO - slow python startup times are a bottle neck

IO currently outputs simple .dat files and immediately calls a python script to plot them. This can lead to many calls to python (once for each dump). The python plotting script is substantially slower than the rest of the code.

Ideally

• output sequentially numbered dat files at each dump
• create a single post-processing script to do the simple visualisations at the end (or, perhaps at the start and have it wait on new dump files) to avoid multiple startup and exits of python
• this also has benefit of paving the way towards keeping dat files for more complicated post-processing / analysis

Transfer notes on each solver to a wiki

as it says in title

CCAPS variants : update the welcome message to indicate which flavour of code used

CCAPS_atm - incorporate new bcs (where appropriate) and drop problem

• incorporate new bc handling where appropriate (hse boundaries NOT necessary?)
• extend to atm specific bc
• add drop problem
• add an outflow for no tangential accn since hse isnt an issue (as oposed to outflow_hse in the other versions of the code)

Does the analytical solution in multigrid test 5 have a typo?

Think its possible (0.001_num / (time+0.001_num)) should be all square rooted, and its just worked out to be fine anyway because you only simulate for small time.

!the analytical solution
  do ix = 1,nx
  do iy = 1,ny
    xx = x(ix) - 0.5_num
    yy = y(iy) - 0.5_num
    r = sqrt(xx**2 + yy**2)
    phi_an(ix,iy) = (phi2-phi1) * (0.001_num / (time+0.001_num)) &
      & * exp(-0.25_num * r**2 / k / (0.001_num + time)) + phi1
  enddo
  enddo

MG / CCAPS - see if OpenMP is worthwhile

e.g. the red-black relax loops in multigrid should be easy enough to put in a parallel do

MG - handle the different levels / grids as a linked list of a derived type of grids