"fix mesh/surface/stress/servo" problems

Submitted by Omid dorostkar on Fri, 03/06/2015 - 10:52

Dear all,

I am about simulation of a box which should be confined vertically , and then sheared horizontally. I have 3 layers, the lower layer (which has 2 sub-layers) is frozen and the upper layer (which has 2 sub-layers as well) is rigid. They are meant to keep the lower layer fixed and move the upper layer (rigid layer) with constant velocity and looking at the micro-mechanical properties of middle layer. For confining vertically , I have used a "fix mesh/surface/stress/servo" to press it down and we know that we cannot use this mesh for move/mesh command to apply shear! Considering these facts, I have got 3 problems now:

1- When I move my stl plane downward and extract the forces on it, it never gets to target value, lot of oscillations, sometimes 0, negative or positive. On the other hand, although I only put force vertically, I see non zero values for horizontal components of forces acting on wall. I tried different "max vel " but still I can not get a permanent stable force on it!

2- The particles go inside each other under pressure! I do not understand this behavior! You can see the attached picture.

2- for next step, when I want to shear my box to see the behavior of middle layer in terms of macro friction coefficient (Fy/Fx, applied to this layer) and due to the problem I discussed in 1st paragraph, I had to make the upper layer rigid, move it horizontally with constant velocity and then calculate forces on the upper wall. It is the simplistic way and not the most precise way, but available one. I am wondering if anybody knows any trick to move a mesh with constant velocity and be able to measure forces on it (force controlled walls). Notice that, this mesh is used in previous stage of experiment (confining).

I have attached a picture, showing my model so far. the two upper layers are rigid and two lower layers are frozen.
Thank you in advance for your responses.

Image icon model configuration206.83 KB
Image icon Particles overlaps are so high!845.69 KB
richti83's picture

richti83 | Fri, 03/06/2015 - 12:42

Hi Omid,

mind sharing your script with us that I can have look ? (you can send it by PM).
I think for shearing the top layer it would be OK to make the top mesh wall frictionless (so only normal forces act) and than use the fix move command to move the particles of the top layer (and not the top wall).

Perform updates of position and velocity for atoms in the group each timestep using the specified settings or formulas, without regard to forces on the atoms. This can be useful for boundary or other atoms, whose movement can influence nearby atoms.

On the other hand, maybe you can use the surface_velocity feature of fix wall mesh or shear feature of primitive wall to not move the triangles (which is not supported for servo walls) to apply a constant surface velocity to the plane. I did not try if this is compatibel with servo wall. If not compatibel I suggest this hack:
1st run until you reached the desired compacting or wallforce,
2nd store position of mesh/surface/stress/servo to variable (varible Pz equal f_[9]The last 3 components output the wall position
3rd make primitive wall at ${Pz} with shear velocity.
4th unfix fix mesh/surface/stress.


Omid dorostkar | Sun, 03/08/2015 - 17:09

I was doing your suggestion like this,

fix cad1 all mesh/surface/stress/servo file pn5.stl type 1 scale 0.001 rotate axis 1. 0. 0. angle 90 move -0.02 0.3123 -0.005 com 0.51 0.3123 0.05 vel_max 0.1 ctrlPV force axis 0. -1. 0. target_val 4.3e+07

fix wallgran all wall/gran model hertz tangential history mesh n_meshes 1 meshes cad1

variable xForce equal f_cad1[1]
variable yForce equal f_cad1[2]
variable zForce equal f_cad1[3]
variable xCom equal f_cad1[7]
variable yCom equal f_cad1[8]
variable zCom equal f_cad1[9]

thermo_style custom step ke c_1 v_xForce v_yForce v_zForce v_xCom v_yCom v_zCom
run 130000 upto

fix omid all wall/gran model hertz tangential history primitive type 1 yplane ${ yCom}

unfix cad1
unfix wallgran
run 30000

but I get this error! "ERROR: Invalid fix ID in variable formula (../variable.cpp:1128)

jessie777 | Fri, 07/19/2019 - 21:37

Dear Omid,

I get the same error here, could you help me with it if you solved?

Thank you!

mschramm | Tue, 07/23/2019 - 18:00

Sorry you are having problems with the servo command.
Could you make a new post topic with the error message you are receiving and the relevant commands you used to get to the error?

Omid dorostkar | Fri, 03/06/2015 - 14:01

Hi Christian,

Thank you very much for your answer. I sent you my script by PM. it is just before shearing.

I will try your suggestion.
Do you have any idea for high overlap of particles. They some times pass through each other!


Omid dorostkar | Fri, 03/06/2015 - 14:09

And one other point christian,

When I move the upper block, the x component of force on the wall changes its sign! it is sometimes +, or - . it is somehow weird! when we move the particles with constant velocity in x direction (+), the x component of force on wall should have constant sign!
it causes many osculations in macro-friction coefficient!


richti83's picture

richti83 | Sat, 03/07/2015 - 13:04

It's in known that rolling friction cdt (type A by Ai et all [1]) is not damped and osscilates until hell freezes over, try epsd2 it's a type C rolling friction modell and damped internally (see [2], Fig 11).

[1] Ai, J.., et.al. : Assessment of rolling resistance models in discrete element simulations
[2] Rahman, M., Katterfeld, A., Schott, D.L., Lodewijks, G.: Influence of the Software on the Calibration process of DEM simulations

Omid dorostkar | Sat, 03/07/2015 - 13:41

I will do that. Actually I am now trying your suggestion for creation of new wall after confinement, to be sheared.
I will be back to discussion with feedback.


Omid dorostkar | Sat, 03/07/2015 - 20:48

I turned out the rolling friction totally, I have now constant y component of force on my servo stl, but when I move rigid body with constant velocity toward +x direction, the x component of force on stl servo fluctuates and changes its sign (+ or -) periodically!

Anyone knows the reason?


JoshuaP | Sun, 03/08/2015 - 11:55

what is your target value? and what you choosed for kp ki kd? can you post the necessary lines and variables for your servo wall?

Omid dorostkar | Sun, 03/08/2015 - 15:30

My target value is 4.3e+07 in cgs units.

Here you can find my command. I did not defined kp ki kd and just used default values!

fix cad1 all mesh/surface/stress/servo file pn5.stl type 1 scale 0.001 rotate axis 1. 0. 0. angle 90 move -0.02 0.3123 -0.005 com 0.51 0.3123 0.05 vel_max 0.1 ctrlPV force axis 0. -1. 0. target_val 4.3e+07

fix wallgran all wall/gran model hertz tangential history mesh n_meshes 1 meshes cad1

Omid dorostkar | Sun, 03/08/2015 - 18:33

Dear Joshuap
Thanks. your suggestion will solve the negative force during moving of driving block or ..... ?

I will try that but I did not get your point.

JoshuaP | Sun, 03/08/2015 - 19:15

your wall is force controlled in y direction so the x-force can be positive or negative and can fluctuate. the y force should be stable and reach the target value more or less.

Omid dorostkar | Sun, 03/08/2015 - 19:39

I can reach to stable target value in y direction perfectly. there is no problem in this regard.

The problem is when I move a bunch of particles which are in a rigid group toward +x direction, the x component of force on the servo is fluctuating extremely, there for for the whole box, I get negative macro friction coefficient! (Fx/Fy of wall).

JoshuaP | Mon, 03/09/2015 - 09:59

as far as I understood you want to see the friction between wall and particles? Maybe try to fix the servo wall after the side pressure is reached, use
fix_modify cad1 integrate stop


Omid dorostkar | Mon, 03/09/2015 - 14:59

Thanks for your comment. I will try that.

deepakpawar.2310 | Sat, 09/07/2019 - 16:17


com in the fix command for servo wall signify the centroid of stl file or centroid of the simulation domain?