Does neumann BC work in AdvectionDiffusionBGKdynamics?

steed188 · October 12, 2020, 2:48am

Hi,
I’m working on thermal convection-diffusion simulation with Boussinesq approximation following the example rayleighBenard3D.

I tried to set BCs like this:

//for velocity, fix velocity = latticeU for inlet, pressure zero-gradient for outlet, bounceback for other walls.
OnLatticeBoundaryCondition3D<T, NSDESCRIPTOR> *nsBC = createInterpBoundaryCondition3D<T, NSDESCRIPTOR>();
nsBC.setVelocityConditionOnBlockBoundaries(nsLattice, inletBox, boundary::dirichlet);
nsBC.addPressureBoundary0P(outletBox, nsLattice, boundary::dirichlet);
setBoundaryVelocity(nsLattice, inletBox, Array<T, 3>(latticeU, 0., 0.));
defineDynamics(nsLattice, roofBox, new plb::BounceBack<T, NSDESCRIPTOR>);
defineDynamics(nsLattice, groundBox, new plb::BounceBack<T, NSDESCRIPTOR>);
defineDynamics(nsLattice, leftsideBox, new plb::BounceBack<T, NSDESCRIPTOR>);
defineDynamics(nsLattice, rightsideBox, new plb::BounceBack<T, NSDESCRIPTOR>);
defineDynamics(nsLattice, frontsideBox, new plb::BounceBack<T, NSDESCRIPTOR>);
defineDynamics(nsLattice, backsideBox, new plb::BounceBack<T, NSDESCRIPTOR>);
setBoundaryVelocity(nsLattice, nsLattice.getBoundingBox(), Array<T, 3>(0, 0, 0));
initializeAtEquilibrium(nsLattice, nsLattice.getBoundingBox(), 1.0, Array<T, 3>(0, 0, 0));
//for temperature, temperature zero-gradient for outlet and frontside and backside, fix temperature for inlet and other walls.
OnLatticeAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR> *adBC = createLocalAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR>();
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, inletBox, boundary::dirichlet );
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, outletBox, boundary::neumann);
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, roofBox, boundary::dirichlet);
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, groundBox, boundary::dirichlet);
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, leftsideBox, boundary::dirichlet);
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, rightsideBox, boundary::dirichlet);	
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, frontsideBox, boundary::neumann);
adBC.setTemperatureConditionOnBlockBoundaries(adLattice, backsideBox, boundary::neumann);

setBoundaryDensity(adLattice, adLattice.getBoundingBox(), 15.0);
initializeAtEquilibrium(adLattice, adLattice.getBoundingBox(), 15.0, Array<T, 3>(0, 0, 0));
setBoundaryDensity(adLattice, leftsideBox, 15.0);	
setBoundaryDensity(adLattice, rightsideBox, 15.0);		
setBoundaryDensity(adLattice, inletBox, 15.0);	
setBoundaryDensity(adLattice, roofBox, 15.0);	
setBoundaryDensity(adLattice, groundBox, 35.0);

I found that the velocity BCs worked well, but the temperature, the Neumann BCs did not work, it was always 1, which is the default value. It was not calculated based on the adjacent grids. Does the Neumann BC work for the advection lattice? Or should I do some specific setting?

By the way, the lattice was setting like:

#define NSDESCRIPTOR descriptors::ForcedMRTD3Q19Descriptor
#define ADDESCRIPTOR descriptors::AdvectionDiffusionD3Q7Descriptor

MultiBlockLattice3D<T, NSDESCRIPTOR> *nsLattice = new MultiBlockLattice3D<T, NSDESCRIPTOR>(nx, ny, nz, new GuoExternalForceSmagorinskyMRTdynamics<T, NSDESCRIPTOR>(parameters.getSolventOmega(), smagoConst));

MultiBlockLattice3D<T, ADDESCRIPTOR> *adLattice = new MultiBlockLattice3D<T, ADDESCRIPTOR>(nx, ny, nz,new AdvectionDiffusionBGKdynamics<T, ADDESCRIPTOR>(parameters.getTemperatureOmega()));

best wishes,

steed188

orestis · October 28, 2020, 8:29am

I think that the

setTemperatureConditionOnBlockBoundaries

does not work with neumann. You can do it manually with the function

addTemperatureBoundary(Box3D domain, BlockLattice3D<T,Descriptor>& lattice,
            boundary::BcType bcType)

on each boundary.

steed188 · October 28, 2020, 9:25am

Dear orestis,
Thank you for your reply. I found that “addTemperatureBoundary” comes from the class
AdvectionDiffusionBoundaryConditionInstantiator3D< T, Descriptor, BoundaryManager >
while not our commonly used
OnLatticeAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR>.
I’ve no idea how to create this kind of boundary.
Is there any example for reference?

best wishes,

steed188

orestis · October 28, 2020, 9:30am

you can use it exactly as IIRC

setTemperatureConditionOnBlockBoundaries

steed188 · October 28, 2020, 10:19am

I’m not sure how should I use it as " setTemperatureConditionOnBlockBoundaries".
I tried like this, but this does not work.

OnLatticeAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR> *adBC = createLocalAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR>();
adBC.addTemperatureBoundary(outlet, adLattice, neumann);

Should I change OnLatticeAdvectionDiffusionBoundaryCondition3D into AdvectionDiffusionBoundaryConditionInstantiator3D? If so, what is BoundaryManager in its template?

best wishes,
steed188

orestis · November 2, 2020, 6:28am

You have different functions to set the boundary conditions, one is

setTemperatureConditionOnBlockBoundaries

which sets BCs on a whole block, but you also have functions like

addTemperatureBoundary2P()

Those add boundary condition cells only on a plane. The 2P here stands for a plane with normal in the “z” direction and pointing in the positive direction “+1”. So the normal is (0, 0, 1). You have the same function with 0P, 0N, 1P, 1N, 2P, 2N. 0 being “x” direction and 1 being “y” direction.
As an example with the local AD BCs you can instantiate a planar BC for a “top” wall like

OnLatticeAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR> *adBC = 
        createLocalAdvectionDiffusionBoundaryCondition3D<T, ADDESCRIPTOR>();

adBC.addTemperatureBoundary2P(top, adLattice, boundary::neumann);

Then you also have similar BCs for edges and corners. setTemperatureConditionOnBlockBoundaries is just a function that calls these functions for you on a complete block.

steed188 · November 2, 2020, 2:52pm

Do you mean addTemperatureBoundaryXX, addTemperatureEdgeXXX, addTemperatureCornerXXX?
I have already tried these manual methods. But these also did not work.
Take the “top” as an example， I tried

adBC.addTemperatureBoundary2P(Box3D(1, nx-2, 1, ny-2, nz-1, nz-1), adLattice, boundary::neumann); //top domain
adBC.addTemperatureEdge1PN(Box3D(0, 0, 1, ny-2, nz-1, nz-1), adLattice, boundary::dirichlet);
adBC.addTemperatureEdge1PP(Box3D(nx-1, nx-1, 1, ny-2, nz-1, nz-1), adLattice, boundary::dirichlet);
adBC.addTemperatureEdge0NP(Box3D(1, nx-2, 0, 0, nz-1, nz-1), adLattice, boundary::dirichlet);
adBC.addTemperatureEdge0PP(Box3D(1, nx-2, ny-1, ny-1, nz-1, nz-1), adLattice, boundary::dirichlet);
adBC.addTemperatureCornerNNP(0, 0, nz-1, adLattice, boundary::dirichlet);
adBC.addTemperatureCornerNPP(0, ny-1, nz-1, adLattice, boundary::dirichlet);
adBC.addTemperatureCornerPNP(nx-1, 0, nz-1, adLattice, boundary::dirichlet);
adBC.addTemperatureCornerPPP(nx-1, ny-1, nz-1, adLattice, boundary::dirichlet);

But during simulation, I found the top domain temperature is always a fixed value. If I have set

setBoundaryDensity(adLattice, top, T0);

Then, top domain will always be T0. If I have not set a temperature for it, then it will always be 1.

By the way, I checked the adjacent nodes (nodes of z= nz-2), which should be the copy source of the neumann BC. And the temperature of the adjacent nodes will gradually approached to T0 or 1 during simulation.
This somehow implies that the temperature of BC nodes is not copied from the adjacent nodes, while it seems that adjacent nodes were calculated from the BC nodes?

with best wishes,
steed188

orestis · November 24, 2020, 7:44pm

Strange. I’ll try to have a look. Would you have a sample of your code somewhere I can read it?

steed188 · November 26, 2020, 8:08am

Here I shared my code.

https://1drv.ms/u/s!Aq8MbuT-tKSEg8Fwzcg7NObYOc36fg?e=l4a2EL

I have deleted all user-defined functions to make the main code clear. The deleted parts are just data processors and does not affect the simulation.

best wishes.

steed188