New Particle Collision Model 0penLB 1.6

[Jijo]

Dear OpenLB community,

I have tried the new particle collision model but I failed in validating the model with experimental studies (Single spherical Particle-wall Collision). The snippet of the code:

// Discretization Settings
int res = 50;
T const charLatticeVelocity = 0.03;

// Time Settings
T const maxPhysT = 0.45; // max. simulation time in s
T const iTwrite = 1E-3; // write out intervall in s

// Domain Settings
T const lengthX = 0.11;
T const lengthY = 0.11;
T const lengthZ = 0.20;

// Fluid Settings
T const physDensity = 1141.;
T const physViscosity = 0.008/1141.;

//Particle Settings
T Radius = 10.0/1000.0;
T const SphereDensity = 2170.;
T const PositionX = 0.055;
T const Positiony = 0.055;
T const PositionZ = 0.16;

Vector<T,3> SpherePosition = {PositionX,PositionX,PositionZ};
Vector<T,3> SphereVelocity = {0.,0.,0.};
Vector<T,3> externalAcceleration = {.0, .0, -9.81 * (1. – physDensity / SphereDensity)};

// Characteristic Quantities
T const charPhysLength = lengthX;
T const charPhysVelocity = 0.66; // Assumed maximal velocity

unsigned contactBoxResolutionPerDirection = 16;
unsigned particleContactMaterial = 0;
unsigned wallContactMaterial = 0;
T youngsModulus = 575e6;
T poissonRatio = 0.25;
T dampingConstant = 0.1;
T coefficientStaticFriction = 0.02;
T coefficientKineticFriction = 0.3;

/// End of the code

The rest is similar to the settlingcube3D example. I do not see any rebound of the particle even if I increase the time step or the grid size. I read the published paper but I could not fix the issue. Any help here?

Regards,

[jan]

Dear Jijo,

as stated in [1], it is necessary to enlarge the particle for the consideration of the contact. You can do this using particle.setField<NUMERICPROPERTIES, ENLARGEMENT_FOR_CONTACT>( /* choose a value */ );. Otherwise it is to be expected that the incompressible fluid slows down the particle so that there is no contact with the wall, especially for high resolutions.

Kind regards,
Jan

[1]: 10.1016/j.partic.2022.12.005

[Jijo]

Hello Jan,

Thank you for your feedback it has been very useful. Just one question I have with the dkt2d example is regarding this part of the code:

//
solidBoundaries.push_back( SolidBoundary<T, DESCRIPTOR::d>(
std::make_unique<IndicInverse<T, DESCRIPTOR::d>>(
cuboid, cuboid.getMin() – 5 * converter.getPhysDeltaX(),
cuboid.getMax() + 5 * converter.getPhysDeltaX()), 2, wallContactMaterial));

//

What I understand is that you gave material number 2 (wall) the properties of wallcontactMaterial but what arguments are you providing here (std::make_unique<IndicInverse<T, DESCRIPTOR::d>>). Also what if I have couple material numbers that I want to consider as walls in that case should I extend the vectorin the same manner?

Thanks,

[jan]

Dear Jijo,

I agree, that part is advanced. Let me break it down:
1. solidBoundaries is an std::vector of solid boundaries (walls). In this case we only have one, but we can have as many as we like. This can be useful if the walls are made of different materials, or if we want to make sure that contacts in a corner or on an edge are detected separately.
2. In the call of the push_back method, we create an instance of a SolidBoundary and add it to the above mentioned vector.
3. The SolidBoundary is created from a unique pointer of an IndicatorF, a latticeMaterial, a contactMaterial and an optional enlargement solely for the contact treatment.

With this in mind, it may be a little easier to understand. The last two entries are the material numbers. The 2 is the one used on the lattice, e.g. in prepareLattice to set the fluid boundary conditions. The second, the wallContactMaterial, is used to identify the wall properties during the contact treatment. As indicated above. Walls can be of different materials but you are likely to use the same boundary condition for the fluid interaction. Therefore the latticeMaterial could be the same but the contactMaterial is different.

The first parameter, which is probably the most complex looking one, is a simple unique pointer to an indicator that describes the geometry of the wall. Here we use the inverse of a cuboid, which describes the fluid domain. This means that the wall contains everything except the fluid. This is done by the first argument of IndicInverse. The other two arguments take the new min and max values. If we didn’t set them manually and just used the min and max of the fluid, the wall would have a thickness of 0, which would probably cause problems with the distance calculations during the contact treatment.

Kind regards,
Jan

[Jijo]

Hello Jan,

Thank you for your comprehensive explanation. I now understand the implementation.I have one final question that I hope you can help me with. I am currently running a convective flow with resolved particles. However, I constantly come across a problem with the particles deforming into the wall although I am using the discrete contact model you have implemented in OpenLB. I have attached a picture for your information:

Before: https://ibb.co/VgdJJ9m
After: https://ibb.co/3B7Yc0K

Code snippet:

particleSystem.defineDynamics<VerletParticleDynamics<T,PARTICLETYPE>>();

// Create solid boundaries for particle interaction
std::vector<SolidBoundary<T, NSDESCRIPTOR::d>> solidBoundaries;

solidBoundaries.push_back( SolidBoundary<T, NSDESCRIPTOR::d>(std::make_unique<IndicInverse<T, NSDESCRIPTOR::d>>(cuboid, cuboid.getMin()-5 * converter.getPhysDeltaX(),cuboid.getMax()+5 * converter.getPhysDeltaX()), 2, wallContactMaterial));

// Create objects for contact treatment
ContactContainer<T, PARTICLECONTACTTYPE, WALLCONTACTTYPE> contactContainer;
// Generate lookup table for contact properties
ContactProperties<T, 1> contactProperties;
contactProperties.set(particleContactMaterial, wallContactMaterial,
evalEffectiveYoungModulus(youngsModulus, youngsModulus,
poissonRatio, poissonRatio),
dampingConstant, coefficientKineticFriction, coefficientStaticFriction);

// Set contact material
for (std::size_t iP = 0; iP < particleSystem.size(); ++iP) {
auto particlesc = particleSystem.get(iP);
setContactMaterial(particlesc, particleContactMaterial);
}

//Check ParticleSystem
particleSystem.checkForErrors();

for ( std::size_t iT = 0; iT < converter.getLatticeTime( maxPhysT ); ++iT ) {
setBoundaryValues(converter, NSlattice,ADlattice, iT, superGeometry);
// Execute particle manager
particleManager.execute<
couple_lattice_to_particles<T,NSDESCRIPTOR,PARTICLETYPE>,
apply_gravity<T,PARTICLETYPE>
>();

// Calculate and apply contact forces
processContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE, ContactProperties<T, 1>>(
particleSystem, solidBoundaries, contactContainer, contactProperties,
superGeometry, contactBoxResolutionPerDirection);

// Solve equations of motion
particleManager.execute<process_dynamics<T,PARTICLETYPE>>();

// Couple particles to lattice (with contact detection)
coupleResolvedParticlesToLattice<T, NSDESCRIPTOR, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(
particleSystem, contactContainer, superGeometry, NSlattice, converter, solidBoundaries);

Any idea why such issue occur? You support is appreciated.
Regards,

[jan]

Dear Jijo,

I’m not sure what I’m seeing in the pictures. Are the particles passing through the wall or are they bouncing back off it?

If they are passing through it:
– Most likely the contact force is too low, perhaps because the Young’s modulus is too low.
– It could also be that the contact handling in your case is not configured correctly. To debug this, I suggest you check if contacts are detected in your case. If they are detected, then the particleContacts or wallContacts of the contactContainer (see https://www.openlb.net/DoxyGen/html/db/d56/structolb_1_1particles_1_1contact_1_1ContactContainer.html) should contain at least one object. You can add this check after coupleResolvedParticlesToLattice.

If they bounce and you’re wondering about the artifact:
If the “automatic” porosity reset doesn’t work for you, you can try to reset the fields manually with resetSuperParticleField (see https://www.openlb.net/DoxyGen/html/da/d81/namespaceolb_1_1particles.html#a9f93e2392819fc95fc4dc52e282893e0) before writing the VTKs.

Best regards,
Jan

[Author]

Posts