Mass conservation in REEF3D_CFD
 

Hello! I'm trying to evaluate if REEF3D_CFD is suitable for some channel flow hydrodynamic problems I need to simulate.

After compiling it, I started by simulating the dam break with obstacle tutorial cases to see how the model run. What I noticed immediately is that in that case the water fraction does not seem to conserve.

At the start of the simulation the water volume is 0.350 m3, and after 2 seconds of simulations, using the default method Level Set method (with 3rd order TVD runge kutta for time discretization, and WENO5 FLUX for convection), it has gone down to 0.267 m3 (24% error).

I tried refining the mesh in half (from 128k cells to 1024k cells), keeping the same Level Set method, but that after 2 sec the fluid fraction still goes down to 0.301 m3 (14% error; 1st order convergence?).

I also tried switching to the VOF formulation (3rd-order TVD Runge-Kutta for time discretization and HRIC modified for convection) in the finer mesh,but the conservation error after 2 sec is 19%.

So, my questions are: where do this conservations errors come from? Are they inherent to the IB method used to represent the obstacle? (probably not, see edit). Any way to reduce them significatively?

Thank you in advance.

Edit:

- I uploaded the test cases to a GitHub repo in case somebody wants to inspect them: https://github.com/nicolasbadano/REEF3D_Tests
- I used the current REEF3D version, cloned from the master branch of the github repo: https://github.com/REEF3D/REEF3D . Seems to be identical to REEF3D 22.04.1
- I did a new test removing the obstacle, to eliminate any influence of the IB method, I still got 14% error after 2 sec
- I attached a table summing up the results for the different tests

Hi Nicolas,

The usual suspect is the grid size, but you seemed to have checked atleast two grids.

Next, you seem to be using obstacles defined by O. Please switch to S. (Will update the manual)


I will look into it.

Thank you Arun, nice to meet you. Thanks for taking the time to look at this.

I did an extra test, concentrating in the case without the obstacle. Hence it's just a simple dam break case, with two discretizations: coarse (Lx/dx = 80; 128k cells) and fine (Lx/dx = 160; 1024k cells).

After 1.5s I'm still getting 8% error on the first mesh and 5% error in the second one). Furthermore, the fine mesh is getting unstable shortly after that time. The conservation errors look too high for such a simple test with such a short simulation time.

I added in attachment a summary for this tests. I also added the cases themselves to the same repository.

Let me know how I can help; I'm genuinely very interested in REEF3D, I'll be happy to contribute in any way I can.

In the meantime I'll repeat the tests with different discretization schemes.

Thanks Nicolas.
Could you please check with k-w (T 30 2) turbulence and D 10 5 for convection as well.
3rd RK should be ok for time.

Hi Nicolas,
I just used your files with the following changes:
control.txt- S in place of O, B 1 0.01

ctrl.txt- T 10 2, T 30 2

I see vol1 reducing from 0.35 to 0.34 from t=0 to t=2s and no further change.

Thank you Arun.

I can confirm that just switching O to S, and the turbulence model to k-omega (T 10 2, T 30 2) greatly reduces the conservation error on the fine mesh I used, from 14% to 4% after 2 sec.

The k-omega model leads to much higher eddy viscosity than the URANS k-e I was using, so the free surface is much smoother (see the attached screenshot). I suspect this is what explains most of the effect.

So, I have three questions:

- Is "T 30 2" really necessary? I can't find any place in the source code (master branch) where it's read. "T 10" is read in "read_control.cpp", line 1474.

- From these experiments, it seems that the free surface has to be "smooth" relative to the cell size for the phase fraction to be conserved accurately. Do you think this is a fair assessment?

- Do you think this limitation could be lifted somehow with some improvement in the model?

Thank you very much again for your time

Hi Nicolas,

maybe have a look at options F 46 and F 47.

By the way, very interesting to see you using VOF. We have not really tested it that much, so good to see you get some results with it. Our VOF algorithm follows a pretty basic approach - as in your case combining HRIC with standard interface compression. I assume there is still a lot of improvement possible if some time and work is invested. Is this research you would be interested in? We also have a PLIC implementation, using the normal calculation based on the level set function. Again, we have not really followed up on this, but I am still intrigued at getting this done and finding out what can be achieved with a properly tested VOF scheme.

You are right, I am not sure why I said that.
I meant to ask you to choose the k-w model, which is done with T 10 2 anyway.

Thank you Hans, this seems to resolve continuity problems perfectly according to my new batch of tests!

I redid the previous simulations, without and with obstacle (with S 10), both with the default options (no Picard iterations) and with Picard iterations, using:

F 46 3
F 47 5 (5 iterations, no idea if this is reasonable or not).

Basically, after 1 sec the simulations without iterations had continuity errors of 4% and 12%. With iterations those errors went down to 0% (according to the reported precision)!

In the attached table I present the summary for this tests. In the attached figure I show results for t=1s; from left to right this are:
- Without obstacle, no Picard Iterations
- Without obstacle, 5 Picard Iterations
- With obstacle, no Picard Iterations
- With obstacle, 5 Picard Iterations

First of all, thank you for developing and letting us use REEF3D. I feel it has much potential.

I use OpenFOAM almost exclusively, but I think that creating conformal meshes is overkill for some consulting open channel flow cases, in which extreme precision is not required. These meshes and solvers have significant costs in terms of man-hours and run time.

Therefore, I've been looking for an alternative model based on the IB method, and REEF3D looks almost perfect for these applications.

I intend to do various validations tests, and this issue with mass conservation was the first thing to grab my attention. I'll report back the results from this tests in case you think they are interesting.

VOF and PLIC look like interesting capabilities for these applications, but I don't think they are absolute requirements. I just tested VOF because I thought that the level set implementation was maybe the source of the errors.

In order to proceed with the validation tests, I still have a couple of issues I need to deal with:

- first and foremost, the manual seems to imply that geometry can be loaded as STL, but I haven't seen a way to do that apart from the 6DoF code. Is this the only intended way?

- I've had many simulations become unstable due to the generation of high velocities near the free surface. I've learnt to deal with this problems OpenFOAM, but I'm not sure how to tackle them yet in REEF3D. Any advice related to this?

Thanks in advance!

There is small inconsistency in how REEF3D handles objects in the domain. In all cases except floating bodies, the objects are defined and entered into the grid from DiveMESH, the mesh generator.
So, to include stl objects into REEF3D, which are not floating bodies, check options S1-S9 in the DiveMESH manual, which you then include in control.txt

This is probably due to overproduction of turbulence at the free surface under complex deformations?
You can activate the free surface turbulence damping algorithm. The paper shows various test cases, but practically:
use T 36 2
set T 37 to 0.07* local water depth (for open channel flow)
For waves, the model is not that sensitive to the values of T 37 and I have tended to use T 37 =0.07 generally or 0.035 for large plunging breaking waves.

Hope that helps :)

Thank you Arun! I'll try these fixes soon.