I wonder if anyone can help? I have extended the sample code segregation2D.cpp in examples\tutorial\shanChenMultiPhase folder to 3D and everything works as expected. I then imported some 3D geometry (as per the permeability tutorial where 1==interface nodes, 0==fluid nodes & 2 ==solid nodes) and defined bounceback on interface nodes, no dynamics on solid nodes and ExternalMomentRegularizedBGKdynamics on fluid nodes. See below:

// Where “geometry” evaluates to 1, use bounce-back.d
defineDynamics(lattice, geometry, new BounceBack<T,DESCRIPTOR>(), 1);
// Where “geometry” evaluates to 2, use no-dynamics (which does nothing).
defineDynamics(lattice, geometry,new NoDynamics<T,DESCRIPTOR>(), 2);
//otherwise on fluid nodes use BGK dynamics
defineDynamics(lattice, geometry, new ExternalMomentRegularizedBGKdynamics<T, DESCRIPTOR>(omega) ,0);

This results in a weird error and I end up with the average density variable being not a number (nan) within the first few iterations. When I take out the bounceback nodes and replace them with no dynamics or BGK dynamics it works as expected. Have I forgot to do something with the bounceback nodes? Does anyone know what the problem may be?

Here’s a possible explanation. Think of fluid nodes which are in contact with bounce-back nodes. In order to compute the Shan/Chen interaction potential, these fluid nodes will access nearest neighbors and compute the density on them. The problem is, density is an undefined value on bounce-back nodes.

In Palabos however, bounce-back nodes behave as if density was well defined, and they simply deliver a density value which has previously been specified by the user. By default, this “virtual density” on Bounce-Back nodes is defined to be 0. If you prefer to have a density of, say, 1 on Bounce-Back nodes, simply create them with a constructor argument as follows:

new BounceBack<T,DESCRIPTOR>(1.)

It’s possible that your instability was triggered by the pressure jump between fluid nodes and solid nodes, and I suggest that you assign to Bounce-Back nodes a virtual density which is of the same order as the initial surrounding fluid density.

Hi, thanks for the advice regarding setting the virtual density to a value equal to the density in the pores. This however does not solve the problem completely. When I run the 3D Shan Chen with porous media it still appears to become unstable. I have left all parameters as they are in the 2D segregation program so am unsure why this does not work. If I change the dynamics from ExternalMomentBGK to regular BGK I do not get the instability but the density field appears to be incorrect i.e. there is no aggregation (segregation) of densities within the pore space as I would expect, instead I get a uniform spreading of the denisity in the pore space within the first 100 time steps.

Is the ExternalMomentBGK required as I am not sure what the external force is in the segregation example? Any other tips or thoughts on the origin of the instability would be greatly appreciated.