How to set the degrees of freedom of the particle in HLBM

[lsflbm]

Hello everyone,
I’m working on an example very similiar with settlingCube3D in HLBM.I put a particle in the microchannel and try to constrain the particle migration along x-axis. Thus, I want to set the degrees of freedom along y and z-axis to zero. Or just simply make the particle moving in a straight line. (particle is allowed freely to rotate).Can someone teach me how to acheive this?
Thank you guys.

Best
Evan

[Nicolas]

Dear Evan,

please always refer to the version of OpenLB you are using, as there are major changes especially in the particle framework.

Assuming, you are using the current one (v 1.5) though, you can do the following with respect to the implementation in settlingCube3D:

The particle motions is taken care of by the

VerletParticleDynamics<T,PARTICLETYPE>

(located in dynamics/particleDynamicsBase.h). This Dynamics-Class is applied to every particle in the example file. Inside, the

velocityVerletIntegration<T,PARTICLETYPE>

(located in functions/particleMotionFunctions.h) is called, which is responsible for modelling Newton’s motion in all cartesion directions.

In order to restrict this calculation to e.g. the z-axis only, just limit the respective position update here to e.g the third vector component (here associated with the z-direction).

You can then incorporate those changes into a new Dynamics Object (similar to the VerletParticleDynamics) and apply this one to your particles in the example file instead.

Best,
Nicolas

[lsflbm]

Dear Nicolas,

Thank you for answering my question!I have already update to v 1.5 and read the new cpp.

I have an idea that I think may work this well. When I want the particle to moving along x-axis,I can reset the particle’s velocity in y and z-axis to 0 before the collide and stream step. Therefore, when the particle is about to move, it only has the velcity in x-axis. Will my idea work? Could the fuction setVelToAnalyticalVel (AnalyticalConst3D< T, T > &) useful to this?

Best,
Evan

[lsflbm]

Dear Nicolas,

Thank you for answering my question!I have already update to v 1.5 and read the new cpp.

I have an idea that I think may work this well. When I want the particle to moving along x-axis,I can reset the particle’s velocity in y and z-axis to 0 before the collide and stream step. Therefore, when the particle is about to move, it only has the velcity in x-axis. Will my idea work? Could the fuction setVelToAnalyticalVel (AnalyticalConst3D< T, T > &) useful to this?

Best,
Evan

[jan]

Dear Evan,

if you want to set the position of the particle, you can access the corresponding field directly. The following should work:

// get particle via id (starts from 0)
Particle<T, PARTICLETYPE> particle = particleSystem.get(particleID);
// set modified position
// you can use a similar method "getField" to get the content of the field
particle.template setField<GENERAL,POSITION>( modifiedPosition );

However, in my opinion it would be best to change the particle’s position before the coupling from the particle to the fluid via couple_particles_to_lattice. Therefore, you might have to separate the execution of the particleManager or even call the specified steps without it.

Best,
Jan

• This reply was modified 2 years ago by jan.

[lsflbm]

Dear Jan,
Thank you for giving suggestion! I tried the method you mention. However, when running the programm, something goes wrong.
I use the settlingCube3d example and add a rectangular poiseuille flow. The particle is a sphere particle. At the firsi time, the programm goes very well. However, at about 5400 step ,the avEnergy and avRho become -nan ,and a few steps later, the position and velocity become very big, like this.
[Timer] step=4800; percent=8; passedTime=306.526; remTime=3525.05; MLUPs=13.9401
[LatticeStatistics] step=4800; t=0.0008; uMax=0.132613; avEnergy=0.000586339; avRho=1.04953
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000326169, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( -2.17393e-12, -2.33405e-12, -1.51876e-11)
[ParticleInfo] |Velocity(m/s)= ( 1.33061, 1.60779e-16, -5.15718e-15)
[ParticleInfo] |Ang. Velocity(°/s)= ( 8.12674e-08, 2.37697e-09, 5.95079e-07)
[ParticleInfo] |Force(N)= ( 4.93697e-09, -6.73433e-18, 6.7951e-20)
[ParticleInfo] |Acceleration(m/s^2)= ( 332.59, -4.53673e-07, 4.57767e-09)
[ParticleInfo] |Ang. acc.(°/s^2)= ( -0.124323, 0.0334074, 6.40006)
[setBoundaryValues] step=4830; maxVel=0.137543
[setBoundaryValues] step=4860; maxVel=0.139692
[setBoundaryValues] step=4890; maxVel=0.141859
[setBoundaryValues] step=4920; maxVel=0.144042
[setBoundaryValues] step=4950; maxVel=0.146242
[Timer] step=4950; percent=8.25; passedTime=314.503; remTime=3497.65; MLUPs=13.9995
[LatticeStatistics] step=4950; t=0.000825; uMax=0.143308; avEnergy=0.000719076; avRho=1.05481
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000360634, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( 7.30257e-13, -2.50328e-12, 1.41442e-11)
[ParticleInfo] |Velocity(m/s)= ( 1.43133, -5.47747e-13, -4.41999e-15)
[ParticleInfo] |Ang. Velocity(°/s)= ( 1.68225e-07, -2.077e-09, 5.91022e-06)
[ParticleInfo] |Force(N)= ( 1.4875e-07, 1.01319e-17, 4.42079e-20)
[ParticleInfo] |Acceleration(m/s^2)= ( 10020.8, 6.82556e-07, 2.97816e-09)
[ParticleInfo] |Ang. acc.(°/s^2)= ( 0.350576, 0.0311973, 27.0941)
[setBoundaryValues] step=4980; maxVel=0.148459
[setBoundaryValues] step=5010; maxVel=0.150692
[setBoundaryValues] step=5040; maxVel=0.152942
[setBoundaryValues] step=5070; maxVel=0.155208
[setBoundaryValues] step=5100; maxVel=0.157491
[Timer] step=5100; percent=8.5; passedTime=322.481; remTime=3471.41; MLUPs=13.996
[LatticeStatistics] step=5100; t=0.00085; uMax=0.154416; avEnergy=0.00087751; avRho=1.06037
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000397766, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( 2.53597e-12, -2.56859e-12, 5.61068e-11)
[ParticleInfo] |Velocity(m/s)= ( 1.54258, -8.64617e-13, -5.10585e-15)
[ParticleInfo] |Ang. Velocity(°/s)= ( -7.01492e-09, -2.54223e-10, 4.43621e-06)
[ParticleInfo] |Force(N)= ( -4.12592e-10, 7.68024e-17, -1.13703e-19)
[ParticleInfo] |Acceleration(m/s^2)= ( -27.7951, 5.17396e-06, -7.65988e-09)
[ParticleInfo] |Ang. acc.(°/s^2)= ( 1.17213, 0.133466, 55.016)
[setBoundaryValues] step=5130; maxVel=0.159789
[setBoundaryValues] step=5160; maxVel=0.162104
[setBoundaryValues] step=5190; maxVel=0.164436
[setBoundaryValues] step=5220; maxVel=0.166783
[setBoundaryValues] step=5250; maxVel=0.169146
[Timer] step=5250; percent=8.75; passedTime=330.475; remTime=3446.38; MLUPs=13.9698
[LatticeStatistics] step=5250; t=0.000875; uMax=19.5993; avEnergy=-nan; avRho=-nan
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000437653, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( -1.33142e-11, -3.37498e-12, 1.37863e-10)
[ParticleInfo] |Velocity(m/s)= ( 1.65302, -1.17893e-12, -1.10465e-14)
[ParticleInfo] |Ang. Velocity(°/s)= ( -7.76043e-07, -1.06017e-07, 0.000170989)
[ParticleInfo] |Force(N)= ( 4.56389e-08, 2.52804e-16, 2.18143e-19)
[ParticleInfo] |Acceleration(m/s^2)= ( 3074.56, 1.70307e-05, 1.46957e-08)
[ParticleInfo] |Ang. acc.(°/s^2)= ( 2.65106, 0.646281, -253.488)
[setBoundaryValues] step=5280; maxVel=0.171525
[setBoundaryValues] step=5310; maxVel=0.17392
[setBoundaryValues] step=5340; maxVel=0.176331
[setBoundaryValues] step=5370; maxVel=0.178758
[setBoundaryValues] step=5400; maxVel=0.1812
[Timer] step=5400; percent=9; passedTime=338.417; remTime=3421.77; MLUPs=14.0595
[LatticeStatistics] step=5400; t=0.0009; uMax=120.978; avEnergy=-nan; avRho=-nan
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000480408, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( -1.44732e-10, -5.16732e-12, 1.76828e-09)
[ParticleInfo] |Velocity(m/s)= ( 1.76688, 6.3939e-11, -1.72505e-14)
[ParticleInfo] |Ang. Velocity(°/s)= ( -1.64761e-05, 6.05737e-08, 0.000779846)
[ParticleInfo] |Force(N)= ( 1.11817e-08, 6.74486e-15, 1.56529e-18)
[ParticleInfo] |Acceleration(m/s^2)= ( 753.279, 0.000454382, 1.05449e-07)
[ParticleInfo] |Ang. acc.(°/s^2)= ( -47.1061, 2.13715, -6333.05)
[setBoundaryValues] step=5430; maxVel=0.183658
[setBoundaryValues] step=5460; maxVel=0.186131
[setBoundaryValues] step=5490; maxVel=0.18862
[setBoundaryValues] step=5520; maxVel=0.191124
[setBoundaryValues] step=5550; maxVel=0.193643
[Timer] step=5550; percent=9.25; passedTime=346.566; remTime=3400.09; MLUPs=13.7041
[LatticeStatistics] step=5550; t=0.000925; uMax=276.399; avEnergy=-nan; avRho=-nan
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.00052615, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( -4.25498e-10, -1.66612e-12, 7.73915e-09)
[ParticleInfo] |Velocity(m/s)= ( 1.89098, 1.05054e-11, -3.86624e-14)
[ParticleInfo] |Ang. Velocity(°/s)= ( -3.72311e-05, 1.88295e-06, -0.0041353)
[ParticleInfo] |Force(N)= ( 7.07657e-08, 1.0805e-14, -3.10701e-18)
[ParticleInfo] |Acceleration(m/s^2)= ( 4767.28, 0.000727903, -2.09311e-07)
[ParticleInfo] |Ang. acc.(°/s^2)= ( -107.366, 2.11181, -15422.1)
[setBoundaryValues] step=5580; maxVel=0.196177
[setBoundaryValues] step=5610; maxVel=0.198726
[setBoundaryValues] step=5640; maxVel=0.201291
[setBoundaryValues] step=5670; maxVel=0.20387
[setBoundaryValues] step=5700; maxVel=0.206464
[Timer] step=5700; percent=9.5; passedTime=354.433; remTime=3376.44; MLUPs=14.1953
[LatticeStatistics] step=5700; t=0.00095; uMax=63.6949; avEnergy=-nan; avRho=-nan
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 0.000574989, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( -1.0073e-09, -1.77059e-10, 5.05713e-08)
[ParticleInfo] |Velocity(m/s)= ( 2.02149, 1.17921e-09, 2.1672e-13)
[ParticleInfo] |Ang. Velocity(°/s)= ( -1.23342e-05, -2.12496e-05, -0.00623013)
[ParticleInfo] |Force(N)= ( 2.03958e-07, -5.41757e-14, -1.14054e-17)
[ParticleInfo] |Acceleration(m/s^2)= ( 13740.1, -0.00364966, -7.68349e-07)
[ParticleInfo] |Ang. acc.(°/s^2)= ( -261.488, -2.67581, 11448.6)
[setBoundaryValues] step=5730; maxVel=0.209073
[setBoundaryValues] step=5760; maxVel=0.211696
[setBoundaryValues] step=5790; maxVel=0.214334
[setBoundaryValues] step=5820; maxVel=0.216986
[setBoundaryValues] step=5850; maxVel=0.219652
[Timer] step=5850; percent=9.75; passedTime=362.293; remTime=3353.53; MLUPs=14.2061
[LatticeStatistics] step=5850; t=0.000975; uMax=23106.1; avEnergy=-nan; avRho=-nan
[ParticleInfo] Particle 0
[ParticleInfo] |Circum radius(m)= 1.58333e-05
[ParticleInfo] |Mass(kg)= 1.4844e-11
[ParticleInfo] |Position(m)= ( 1.09741e+41, 2.5e-05, 5e-05)
[ParticleInfo] |Angle(°)= ( 26.2179, -14.427, 50.8281)
[ParticleInfo] |Velocity(m/s)= ( 4.11529e+46, -1.39122e+41, -9.36546e+38)
[ParticleInfo] |Ang. Velocity(°/s)= ( 3.13484e+42, -7.97655e+45, 1.18674e+48)
[ParticleInfo] |Force(N)= ( 0, 0, 0)
[ParticleInfo] |Acceleration(m/s^2)= ( 0, 0, 0)
[ParticleInfo] |Ang. acc.(°/s^2)= ( 0, 0, 0)
[setBoundaryValues] step=5880; maxVel=0.222333
[setBoundaryValues] step=5910; maxVel=0.225028
[setBoundaryValues] step=5940; maxVel=0.227737
[setBoundaryValues] step=5970; maxVel=0.23046
[setBoundaryValues] step=6000; maxVel=0.233196

Could you plz tell me what’s going on?

Best
Evan

[jan]

Dear Evan,

due to the given output I assume that the lattice velocity () becomes too large. Try a smaller velocity (1/10 or even less just for testing). If the simulation is stable with the smaller velocity, then you probably need to change your discretization (Δt and Δx) in order to make the simulation work with the original velocity.

Also, do you still apply the gravitational force (it’s called external acceleration in the example)? This could cause an erroneous particle velocity.

Best,
Jan

[Nicolas]

Dear Evan,

before considering any particle related issue, please ensure that the described problem does not appear in the particle-free flow after the same duration already. If it does, you need to respect stability and accuracy limits first!

Best,

Nicolas

[lsflbm]

Dear Jan,
I tried to set the lattice velocity to 1/10 but it didn’t work, it still goes quite large at about 10% of the whole time as before. Besides, I didn’t add gravity in this model. It’s just a sphere particle in a rectangular channel with rectangular poiseuille flow.
In addition, I have a lilttle problem of what you mentioned that change the particle’s position before the coupling from the particle to the fluid via couple_particles_to_lattice. Did you mean change the position after couple_lattice_to_particles and before the couple_particles_to_lattice ? I seperate the execute to couple lattice, change postion, couple patcicle. But it still couldn’t work. Could you tell me more details?

Thank you again for answering me!

Best,
Evan

[lsflbm]

Dear Nicolas,

Thank you for answering my question.

However, I just start leaning LBM and Openlb for not a long time and I couldn’t really get the point of what you said. Could you plz tell me more details? Do you mean that the rectangular poiseuille flow is a particle-free flow? What else could I do besides confirming the CFL condition?

Thank you!

Best,
Evan

[mathias]

Dear Evan,

please consider to come to our next spring school. In that framework I am sure we can answer many of your questions.

Best
Mathias

[Author]

Posts