[jan]

Dear avrachan1,

what do you mean by “I was representing the arbitrary shaped particle by decomposing it into smaller cuboids. This method doesn’t scale well.”? If you use the function addResolvedArbitraryShape3D then it shouldn’t add any additional costs during the simulation, because that leads to the use of a cached signed distance in the background. The interpolation using that cache doesn’t depend on the original geometry.

For complex shapes, I usually use the combination of geometries using signed distance functions or STLs. Depending on what shapes you need, you could try to use a superellipsoid or try to introduce a new indicator that serves your needs, but the main options that are currently implemented are the STLreader or the combination of primitives.

Best regards,
Jan

[jan]

Dear jakschee,

BLOCKFUNCTOR is the functor on the actual core structures. Therefore, for the SuperLatticeMomentumExchangeForce the BlockLatticeMomentumExchangeForce would be used on each separate block. Please refer to the user guide and checkout the difference between super and block lattice.

Periodic is a vector of bools, that identifies which directions are periodic or not. iP0 allows to set an offset at which the particles loops should start, I suggest using the default value (0) here. ignoredMaterials allows you to set lattice materials that are ignored during the evaluation of the functor, meaning that cells having that material number are not considered. f refers to a function that can be passed in case some modification or output is necessary. All of these have default values, so you actually don’t have to define them, if you don’t need that functionality.

Best regards,
Jan

[jan]

Dear araeli,

Do I understand that your question is about the layer of zeros around the actual geometry you imported?

– Is this behavior normal?
Yes, this is normal. Without periodic boundaries, OpenLB sets a layer of 0 on the outside. (Feel free to check some examples).

– Is the excess always just one cell?
Yes, in many examples an extra layer is added to the geometry using IndicatorLayer3D (e.g. cylinder3d, line 336 in cylinder3d.h).

– Shift:
I think this is caused by the discretization? In OpenLB you only have the information if a certain point is inside or outside, and this information seems to match your geometry, doesn’t it?

– Is there a way to force the boundary mesh to adhere to the vertices of the porous medium once they are given as voxels, thus theoretically already cubic and easily adaptable in structure?
Is this related to the shift above? As mentioned above, this seems to be caused by the discretization. Maybe you can make the simulation domain “better” match your geometry (does the STL start at the origin or is it shifted a bit?). However, this discretization error should get better with increasing resolution. If it is your intention to keep a structure made up of smaller cubes, you could also consider using many IndicatorCuboid3D to create in the OpenLB case, like this it might be easier to align them with the lattice.

Best regards,
Jan

[jan]

Hi all,

@nipinl thanks again for the valuable input.

Just a few additional thoughts on this (or maybe just reiterating what nipinl said in other words):

If you don’t want the periodic boundary, then you should turn it off, as nipinl mentioned, by setting everything to false cuboidGeometry.setPeriodicity(false, false, false); or just remove the line (it’s off by default).

I haven’t tested it, maybe you’re lucky and the current implementation already sets the walls correctly. However, you should make sure it does by checking the material numbers (the wall should be 2). If it doesn’t, then you should add a single layer of material number 2 to all 4 faces in prepareGeometry, similar to the inflow and outflow. Note that the superGeometry.clean(); has the potential to mess up the material numbers, so if you need to change something, you might want to turn it off first. Otherwise the changes are sometimes harder to evaluate with clean enabled.

Best regards,
Jan

[jan]

Hello all,

if you believe that bounce back might be the cause, then I’d suggest to try another boundary condition, see examples/laminar/poiseuille2d for a few options and to compare them.

@nipinl it’s great that you already considered the stable limits for the parameters. When you consult the paraview output, can you identify where the problem originates? For example some corner cell or at the inlet or outlet? Another example: If you have the inlet or outlet right next to a “wall” cell, then perhaps this causes problems.

In general I suggest the following:
– Where is the origin? Is there some pulsing? (then consider increasing the startup time or the simulated physical time frame) Is it located at a certain location in the geometry? Etc.
– Check the geometry, are the material numbers as you expect them? Are there any areas with potential conflict (like the example I gave above)?

Once we know a little more details, perhaps we can find the underlying problem easier.

Best regards,
Jan

[jan]

Dear mdiazv64,

You should check that the material numbers are correct (that the inlet, outlet, and walls are where you want them to be). Also, as nipinl pointed out, you should try to reduce the pressure drop. Maybe it’s too high.

@nipinl: Thanks for the valuable suggestion.

Best regards,
Jan

[jan]

Dear Yueyq,

you can try use particles::access::setInactive(particle);, when the particles touch the wall (or any other condition). However, then they’re not going to move anymore.

> Do you think it’s feasible if I set particles velocity and acceleration to be 0, and fixed their positions as long as particles contact the wall?
If you set the velocity and acceleration to 0 and keep the position fixed, then the particles should stay at the wall. Is that what you want? That the particles touch the wall and stick there in any case?

Best regards,
Jan

[jan]

Dear yueyq,

the example only works in parallel mode. The particle system currently doesn’t run on GPUs. However, you can run it in parallel using multiple CPUs. Please refer to the user guide (https://www.openlb.net/wp-content/uploads/2024/06/olb_ug-1.7r0.pdf), see “10.9. Lesson 9: Run Your Programs on a Parallel Machine”.

Best regards,
Jan

Edit:
Yes, the logic of VerletParticleDynamicsWallCapture is only for subgrid scale particles. In my opinion, you’d need to somehow model cohesive forces on the resolved surface for surface resolved particles.

• This reply was modified 10 months ago by jan.

[jan]

Dear araeli,

I believe the following examples may help you:
– examples/laminar/cylinder3d
– examples/particles/bifurcation3d
– examples/adsorption/microMixer3d
– examples/turbulence/aorta3d
– examples/porousMedia/city3d

In general, it is also helpful to check the documentation:
– User guide: https://www.openlb.net/wp-content/uploads/2024/06/olb_ug-1.7r0.pdf (3.4. Reading STL-files)
– Doxygen: https://www.openlb.net/DoxyGen/html/d2/d60/classolb_1_1STLreader.html#a0c802e7d7c251615de04b2610697d0d5

The voxelSize should usually be converter.getConversionFactorLength(). However, the stlSize (3. parameter) depends on your STL. If your STL is in meter, it should be 1. If your STL is in millimeter, it should be 1e-3, etc. Therefore, 1e-6 should be correct in your case.

What error do you get exactly and what are the dimensions of your geometry? Are you certain that each dimension is resolved sufficiently?

Best regards,
Jan

[jan]

Dear Bobbie,

if your 2D structure can be converted to a VTI, then you can use an approach such as this, which uses the IndicatorBlockData2D:

    std::shared_ptr<BlockVTIreader3D<T, T>> vtiReader(
        new BlockVTIreader3D<T, T>(vtiFile, "Elevation"));
    vtiReader->printInfo();
    auto           cuboidSample     = vtiReader->getCuboid();
    T              deltaRsample     = cuboidSample.getDeltaR() * sourceScale;
    Vector<int, 3> extentSample     = cuboidSample.getExtent();
    Vector<T, 3>   originSamplePhys = cuboidSample.getOrigin() * sourceScale;
    Vector<T, 3>   extentSamplePhys = {deltaRsample * T(extentSample[0] + 0.5),
                                       deltaRsample * T(extentSample[1] + 0.5),
                                       deltaRsample * T(extentSample[2] + 0.5)};
    std::shared_ptr<IndicatorBlockData2D<T>> ind(
        new IndicatorBlockData2D<T>(vtiReader->getBlockData(), extentSamplePhys,
                                    originSamplePhys, deltaRsample, false));

Best regards,
Jan

[jan]

Dear nipinl,

did you run the case from the example or another one? For me, it works as described there, even without a cluster. You should try to run that first and get accustomed to the parameters of the setup using that example.

The example case also comes with an explanation of the command line arguments:
<filename> <arrayname> <scaling-factor> <time-scaling-factor> <resolution> <pressure_drop>
So in the example you mentioned, rock.vti would be the file with the geometry, "Tiff Scalars" is the name of the array with the wanted data in that file.
Since the files are usually in different units, a scaling factor is provided. In this case the scaling factor is 2.5e-6, but you can use anything else. A pretty side effect is that you can actually artificially increase and decrease your pore size. A similar scaling factor is given for time. In the example 1 is suggested, but this can be changed if you want to have a longer or shorter simulation. 200 refers to the resolution (it resolves the y-length (extent[1]) with 200 cells -> physDeltaX is calculated from it). Finally, you have to option to define a pressure drop, which is set to 1.

What do you mean by “I replace 200 with number of slices”, how many slices do you have? It sounds like the resolution you’re using is much too high. Just start with a resolution of 200 or perhaps even 100 first and then perform grid independence studies to find the necessary resolution.

Best regards,
Jan

[jan]

Dear yueyq,

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));

the above part takes the inverse of the whole computation domain, which means each side is regarded as a wall. If you solely want the bottom to act as a wall, then I suggest that you introduce a IndicatorCuboid instead and place it at the bottom wall. Be sure to place it exactly where the bottom wall is and give it a sufficient thickness, because the particles will overlap it.

That’s also part of the reason why the SolidBoundary in the example is extended. If you’d have a wall with a infinitesimal thickness, it would lead to problems during the contact treatment, as we want to calculate the overlap volume.

Best regards,
Jan

• This reply was modified 10 months, 1 week ago by jan.

[jan]

Dear Rookie,

I’m not familiar with the legacy framework, but if I recall correctly, the problem you’re experiencing occurs when a large number of particles are used, or actually when a large number of particles cross a boundary to another block. I believe the new particle framework handles these cases better, but I don’t know of a quick fix for the legacy framework.

Best regards,
Jan

[jan]

Dear Rookie,

does that mean that removing p->setCuboid(newCuboid) seems to fix the simulation? Could you explain to me why “particle allocation was no longer necessary”? When is the particle then “transferred” to the neighboring block?

Best regards,
Jan

[jan]

Dear yueyq,

Without having the specific details of your code, I can only make assumptions about the issues you’re facing. Here are a couple of points to consider:

1. It’s possible that you’ve defined the solidBoundaries incorrectly. Ensure that the solidBoundaries are not set for the entire surrounding part, but only for the single wall at the bottom of the container. Updating this might resolve the issue with forces being exerted by the left and right boundaries.

2. Double-check how you’ve set the periodicity. Also, incorrect definition of the walls might cause the force from the wall contact to override the forces from the fluid. This could be the reason why there’s no momentum exchange between the particle and the fluid after the particles leave and re-enter the domain.

I hope these suggestions help you resolve the issues in your code. If you need further assistance or have any other questions, please don’t hesitate to ask.

Best regards,
Jan

[Author]

Posts