Validation of HLBM

[avrachan]

Dear Community,

I have been trying to validate the HLBM module of OpenLB with analytical result of the low reynolds number flow through a periodic array of spheres, measuring the friction coefficient (https://ibb.co/6QvQg3P). [1]

I am unable to match the values as shown in the linked figure (https://ibb.co/zVDpTLW2), the friction coefficient is consistently underestimated.

It would be great if someone can point at some work where such a validation in done using OpenLB or a partially saturated method as implemented in OpenLB.

I have also linked the code here.

Thanks,

https://gist.github.com/vavachan/f1e3ca2d8b8a035f4e1dc7fbd2a590d6

[1] Ladd, A. J. C. Numerical Simulations of Particulate Suspensions via a Discretized Boltzmann Equation. Part 2. Numerical Results. J. Fluid Mech. 1994, 271, 311–339. https://doi.org/10.1017/S0022112094001783.

• This topic was modified 1 month, 2 weeks ago by avrachan.

[stephan]

Dear avrachan,

thank you for your post.

– Could you please explain how the quantities in the linked plots are defined exactly?
– Which range of porosity / permeability would you like to model?
– What are the mesh sizes you consider?
– How is the set up of the domain? How large is it? How many spheres are in it? Is there an inflow and an outflow region? Which boundary conditions are you using? What is Reynolds, Mach, and relaxation time exactly?

In this dissertation (https://dx.doi.org/10.5445/IR/1000161726), Section “4.5 Homogenized Navier–Stokes equations” describes a test case based on laminar flow through an array of spheres.

BR
Stephan

BR
Stephan

[avrachan]

Hi Stephan,

I define an L X L X L cube with periodic BC in all direction. I place a sphere at the center whose radius is varied to change the packing fraction of the system (vol-sphere/L^3).

All fluid nodes are subjected to a arbitrary small force in the y direction. I have used forcing::ForcedHLBM and defined the FORCE field in the lattice.

Upon reaching steady state, the force acting on the sphere is equal to the total force applied to the fluid nodes.

Here we measure the force acting on the sphere in the y direction (F), the mean velocity in the whole system in y direction (U).

Then I plot F/U*(1/6 pi eta a) where eta is the shear viscosity and a is the radius of the sphere, w.r.t to the packing fraction.

I have tried L = 16,32,64,128. I tried packing fractions from ~0 to 0.4. I have tried the relaxation tau close to 0.5 to 1.0. The forcing magnitude from 1e-5-1e-7 ( Re < 1).

The steady state values of the friction coefficient is largely independent of the choices of tau, forcing magnitude and system size (as it should be).

I feel like I am missing something!

Thanks for reading.

[avrachan]

dx=1 and dt=1

Changing L while keeping the packing fraction same changes the resolution, so I changed different parameters but the results are consistently different from what I expect.

I am reading the thesis, and will see how I can compare with it.

Thanks.

[avrachan]

Hi Stephan,

I have attached here the information about trying the drag on one packing fraction at 0.4.

The correct value of the friction coefficient should be around 20 according to Ladd et al [1], but I am getting a value of ~14.

I have realized that I have to reduce the value of tau to near 0.5 to get a value closer to what the Ladd et al shows.

Since this is a fundamental test, I am unable to resolve the issue.

`[UnitConverter] —————– UnitConverter information —————–
[UnitConverter] — Parameters:
[UnitConverter] Resolution: N= 15
[UnitConverter] Lattice velocity: latticeU= 1.64572e-07
[UnitConverter] Lattice relaxation frequency: omega= 1.96078
[UnitConverter] Lattice relaxation time: tau= 0.51
[UnitConverter] Characteristical length(m): charL= 14.6265
[UnitConverter] Characteristical speed(m/s): charU= 1.64572e-07
[UnitConverter] Phys. kinematic viscosity(m^2/s): charNu= 0.00333333
[UnitConverter] Phys. density(kg/m^d): charRho= 1
[UnitConverter] Characteristical pressure(N/m^2): charPressure= 0
[UnitConverter] Mach number: machNumber= 2.85047e-07
[UnitConverter] Reynolds number: reynoldsNumber= 0.000722134
[UnitConverter] Knudsen number: knudsenNumber= 0.000394729
[UnitConverter] Characteristical CFL number: charCFLnumber= 1.64572e-07
[UnitConverter]
[UnitConverter] — Conversion factors:
[UnitConverter] Voxel length(m): physDeltaX= 1
[UnitConverter] Time step(s): physDeltaT= 1
[UnitConverter] Velocity factor(m/s): physVelocity= 1
[UnitConverter] Density factor(kg/m^3): physDensity= 1
[UnitConverter] Mass factor(kg): physMass= 1
[UnitConverter] Viscosity factor(m^2/s): physViscosity= 1
[UnitConverter] Force factor(N): physForce= 1
[UnitConverter] Pressure factor(N/m^2): physPressure= 1
[UnitConverter] ————————————————————-`

[1] Ladd, A. J. C. Numerical Simulations of Particulate Suspensions via a Discretized Boltzmann Equation. Part 2. Numerical Results. J. Fluid Mech. 1994, 271, 311–339. https://doi.org/10.1017/S0022112094001783.

[stephan]

Dear avrachan,

can we fix a single set of physical parameters you are trying to simulate?
Like, a set of fixed Re, volume fraction, L.
Then choose some stable set of discrete parameters (for example with the help of Krueger et al. book, there’s a section on choosing parameters).
Then simulate that for a range of resolutions in diffusive scaling (Ma->0, where dt scales as dx^2).

Could you please try that? Does it work then?

If that approximates a single configuration in your target plot, you can move on to the next one.

BR
Stephan

[avrachan]

Hi Stephan,

Thank you for your reply.

I have tried various combination of parameters.

For example, I have fixed the packing fraction to be 0.3, and applied a very small forcing field ~ 10^{-10} and varied the relaxation time tau.

I find a strong dependence on the value of the drag coefficient on the value of tau. F/v * [6 pi et a]^(-1), is supposed to be independent of tau. This has been reported in [1] section 3, where they looked at a packing fraction below 0.1. [1] is the only comparative study of partially saturated method (as implemented in OpenLB) with Momentum exchange method as discussed in Ladd [2] where drag force on periodic array of spheres in calculated. At the packing fraction they have considered, the calculated value from OpenLB agrees with analytical calculation. The disagreement between analytical values and values calculated from my code increases at higher packing fractions, where I have found that I have to chose very small tau (close to 0.5) to get values close to the analytical one.

As shown in the previous post, the Re, Ma and Kn are all very small. I find the strongest dependence on Kn, due to dependence on tau.

Analytical value of the drag coefficient have been reproduced by [3] (Table 1) for the TRT dynamics and momentum exchange method.

This makes me wonder, if the issue is related to the HLBM and the BGK dynamics.

PS: I have linked the code in the very first post.

[1] Rettinger, C.; Rüde, U. A Comparative Study of Fluid-Particle Coupling Methods for Fully Resolved Lattice Boltzmann Simulations. Computers & Fluids 2017, 154, 74–89. https://doi.org/10.1016/j.compfluid.2017.05.033.

[2] Ladd, A. J. C. Numerical Simulations of Particulate Suspensions via a Discretized Boltzmann Equation. Part 2. Numerical Results. J. Fluid Mech. 1994, 271, 311–339. https://doi.org/10.1017/S0022112094001783.

[3] Bogner, S.; Mohanty, S.; Rüde, U. Drag Correlation for Dilute and Moderately Dense Fluid-Particle Systems Using the Lattice Boltzmann Method. International Journal of Multiphase Flow 2015, 68, 71–79. https://doi.org/10.1016/j.ijmultiphaseflow.2014.10.001.

• This reply was modified 3 weeks, 6 days ago by avrachan.

[stephan]

Dear avrachan,

as far as I know, there is no issue with HLBM and the BGK dynamics are validated.
The problems you are facing are very likely due to the choice of parameters.

1) I would suggest to paste your code here, if you like.
2) Please try to reproduce only one configuration from one reference at a time.
3) Then follow the book of Krüger et al. to choose your parameters. I would propose to choose a real physical characteristic length (the side length of your problem domain for example). Then choose the grid spacing and the viscosity and the real physical flow speed (or the Reynolds number alternatively). Then set the lattice velocity to be smaller then 0.1. This will give you the relaxation time and the time step size delta t.
4) Double-check the settings and the results with respect to this one particular physical setting until the solver converges.
5) Please print the converter of this particular setting again.

BR
Stephan

[Author]

Posts