Capillary pressure between two particles

Hello Seungcheol,

after executing the script in command line, you are using play button in
GUI controller? To track particle displacement, add

from yade import plot

in the beginning of your script and add PyRunner engine at the bottom of
engines:

O.engines = [
   ...
   PyRunner(iterPeriod=5, command="plot.addData(
b0displ=O.bodies[0].state.displ(), b1displ=O.bodies[1].state.displ() )" )
]

change 5 to suitable value. See [2] for docs and [1] why 'plot' is used to
store values. If you add more values, it is maybe better to all the
plot.addData(..) write in separate function and call only that one in
PyRunner, see [2].

Another possibility how to see small movements of particles is to
exaggerate the displacement in view window. In the controller window, go to
Display tab and change dispScale values to something > 1.

I don't use capillary laws, so I can't help you with the other question

HTH
Jan

[1] https://yade-dem.org/doc/user.html#running-python-code
[2] https://yade-dem.org/doc/yade.plot.html

2013/7/24 Seungcheol Yeom <email address hidden>

> New question #232941 on Yade:
> https://answers.launchpad.net/yade/+question/232941
>
> Hello all,
>
> I am trying to run a simulation to see whether two particles are pulling
> each other due to the capillary pressure.
> However, it seems that two particles are not moving at all graphically.
> Can anyone tell me how to track either a capillary pressure or a
> displacement in order to check they are moving?
>
> Here is the script I am using:
>
> from yade import qt
>
> r = 1e-4 #particle radius
> h = 1e-5 #praticle distance
>
> #create two sphere paticles#
> O.bodies.append([
> utils.sphere(center=(0,0,0),radius=r,fixed=False),
> utils.sphere((0,0,2*r+h),r)
> ])
>
> #define engines#
> O.engines=[
> ForceResetter(),
> InsertionSortCollider([Bo1_Sphere_Aabb()]),
> InteractionLoop(
> [Ig2_Sphere_Sphere_L3Geom()],
> [Ip2_FrictMat_FrictMat_CapillaryPhys()],
> [Law2_ScGeom_FrictPhys_CundallStrack(neverErase=True)]
> ),
>
> Law2_ScGeom_CapillaryPhys_Capillarity(capillaryPressure=10000,label='Cap'),
> NewtonIntegrator(damping=0.4,gravity=(0,0,0))
> ]
>
> Cap.createDistantMeniscii=True
> O.run(1,True)
> Cap.createDistantMeniscii=False
>
> O.dt=0.5*PWaveTimeStep()
> qt.View()
> O.saveTmp()
>
> FYI, I am using yade-daily for this script
> I sincerely would appreciate for your help.
>
> Seungcheol
>
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

You can also investigate the script I wrote for getting a force-distance-plot.
There data is stored in python lists and plotted with pylab/matplotlib.
If you want, you can modify this script and add the displacement O.bodies[X].state.displ(), like Jan suggested.

hope it helps,

christian

#!/usr/bin/python
# -*- coding: utf-8 -*-

### properties:
density = 2650
shear_modulus = 2e5
poisson_ratio = 0.15
young_modulus = 2*shear_modulus*(1+poisson_ratio)
suction = 3000
velocity = 0.1 #in x-direction
a=True

### define a material:
mat = FrictMat(young=young_modulus,poisson=poisson_ratio,density=density,frictionAngle=1)
O.materials.append(mat)

### append 2 spheres:
R=.0003
d=.000015 #positive for overlaps, negative for distant
O.bodies.append(sphere([0,0,0],radius=R, material = mat, fixed=True))#id 0
O.bodies.append(sphere([2*R-d,0,0],radius=R, material = mat, fixed=True))#id 1

### define engines:
O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_MindlinCapillaryPhys(betan=0,betas=0,label='ContactModel')],
  [Law2_ScGeom_MindlinPhys_Mindlin(neverErase=True,label='Mindlin')]
  ),
 Law2_ScGeom_CapillaryPhys_Capillarity(capillaryPressure=suction,label='CapPhys'),
 NewtonIntegrator(damping=0,label='integrator'),
 PyRunner(iterPeriod=1,command='a=getData(a)')
]
#CapPhys.dead=True
O.bodies[1].state.vel=Vector3(velocity,0,0)
O.dt=1e-6

capForce=[]
normForce=[]
normCapForce=[]
distance=[]

def getData(a):
 distance.append(O.bodies[1].state.pos[0]-O.bodies[0].state.pos[0]-2*O.bodies[0].shape.radius)
 try:
  i = O.interactions[0,1]
  if i.isReal:
   capForce.append(-1*i.phys.fCap[0])
   normForce.append(i.phys.normalForce[0])
   normCapForce.append(i.phys.normalForce[0]-i.phys.fCap[0])
  else:
   if a:
    print 'interaction lost at step',O.iter
    a=False
   capForce.append(0)
   normForce.append(0)
   normCapForce.append(0)
 except IndexError:
  capForce.append(0)
  normForce.append(0)
  normCapForce.append(0)
 return a

if d <= 0.0:
 CapPhys.createDistantMeniscii=True
 O.run(1,True)
 CapPhys.createDistantMeniscii=False

O.run(600,True)

### make nice figure:

from pylab import * #overwrites angle, box, ... see 0-generate.py
rc('text',usetex=True)
rc('font',**{'family':'serif','serif':['Computer Modern']})

from matplotlib.font_manager import FontProperties
font = FontProperties()
font.set_size('x-large')

figure()
axhline(color='gray')
axvline(color='gray')
plot(distance, normCapForce, label = 'normal + capillary force')
plot(distance, normForce, 'g--', label = 'normal force')
plot(distance, capForce, 'r--', label = 'capillary force')
xlabel('distance in [$m$]',fontproperties=font)
ylabel('force in [$N$]',fontproperties=font)
ylim([-0.00015,0.00025])
xlim([-0.00002,0.00005])
#grid()
legend()
#savefig('force-dist-plot-cap+norm.png')
show()

You can also use displacement scaling to see small displacements better - https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.OpenGLRenderer.dispScale .

Hello all,

First of all, thanks for all your response. I have learned a lot from you.
By the way, I have modified the script as shown below:

 from yade import plot,qt

r = 1e-4 #particle radius
h = 1e-5 #praticle distance

#create two sphere paticles#
O.bodies.append([
   utils.sphere(center=(0,0,0),radius=r,fixed=False),
   utils.sphere((0,0,2*r+h),r)
])

#define engines#
O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_CapillaryPhys()],
      [Law2_ScGeom_FrictPhys_CundallStrack(neverErase=True)]
   ),
   Law2_ScGeom_CapillaryPhys_Capillarity(capillaryPressure=10000,label='Cap'),
   NewtonIntegrator(damping=0.4,gravity=(0,0,0)),
   PyRunner(iterPeriod=100,command='addPlotData()')
]

Cap.createDistantMeniscii=True
O.run(1,True)
Cap.createDistantMeniscii=False

def addPlotData():
    plot.addData(i=O.iter,b0displ=O.bodies[0].state.displ(),b1displ=O.bodies[1].state.displ())

plot.plots={'i':('b0displ')}
plot.plot()
plot.live=True
plot.autozoom=True
O.dt=0.5*PWaveTimeStep()
qt.View()
O.saveTmp()

When I ran it, I am not able to see any displacements based on the plot. I hoped I can see the displacement as a function of the number of iteration but it seems something went wrong.
Am I making something wrong?

Thanks for your help!

Seungcheol

Hello Seungcheol,

probably just a very stupid question, but do you actually run the
simulation (I already asked last time, if you run it in GUI by play
button)? if yes, how your stored data looks like? you can access it in
ipython by
plot.data
command
cheers
Jan

2013/7/24 Seungcheol Yeom <email address hidden>

> Question #232941 on Yade changed:
> https://answers.launchpad.net/yade/+question/232941
>
> Seungcheol Yeom posted a new comment:
> Hello all,
>
> First of all, thanks for all your response. I have learned a lot from you.
> By the way, I have modified the script as shown below:
>
> from yade import plot,qt
>
> r = 1e-4 #particle radius
> h = 1e-5 #praticle distance
>
> #create two sphere paticles#
> O.bodies.append([
> utils.sphere(center=(0,0,0),radius=r,fixed=False),
> utils.sphere((0,0,2*r+h),r)
> ])
>
> #define engines#
> O.engines=[
> ForceResetter(),
> InsertionSortCollider([Bo1_Sphere_Aabb()]),
> InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom()],
> [Ip2_FrictMat_FrictMat_CapillaryPhys()],
> [Law2_ScGeom_FrictPhys_CundallStrack(neverErase=True)]
> ),
>
> Law2_ScGeom_CapillaryPhys_Capillarity(capillaryPressure=10000,label='Cap'),
> NewtonIntegrator(damping=0.4,gravity=(0,0,0)),
> PyRunner(iterPeriod=100,command='addPlotData()')
> ]
>
> Cap.createDistantMeniscii=True
> O.run(1,True)
> Cap.createDistantMeniscii=False
>
> def addPlotData():
>
> plot.addData(i=O.iter,b0displ=O.bodies[0].state.displ(),b1displ=O.bodies[1].state.displ())
>
> plot.plots={'i':('b0displ')}
> plot.plot()
> plot.live=True
> plot.autozoom=True
> O.dt=0.5*PWaveTimeStep()
> qt.View()
> O.saveTmp()
>
> When I ran it, I am not able to see any displacements based on the plot. I
> hoped I can see the displacement as a function of the number of iteration
> but it seems something went wrong.
> Am I making something wrong?
>
> Thanks for your help!
>
> Seungcheol
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Hello Jan,

I ran the simulation by hitting a play button in GUI.
By typing in plot.data in the command, I am seeing the numbers only.
It seems like it is the iterperiod because the number is increased by 100.
Am I missing something?
Thanks a lot Jan!

Seungcheol

Hello Seungcheol,

plot.data is python dictionary. Keys defined by plot.addData function and
values are lists of numbers. So improved suggestion is to type
plot.data['b0displ']
or
plot.data['b1displ']

If it is full of zeros, then you have really no displacement, if the
numbers are nonzero, then you have some displacement. After this it should
be clear where the problem is (only displaying for example or something
inside the simulation)
cheers
Jan

2013/7/24 Seungcheol Yeom <email address hidden>

> Question #232941 on Yade changed:
> https://answers.launchpad.net/yade/+question/232941
>
> Seungcheol Yeom posted a new comment:
> Hello Jan,
>
> I ran the simulation by hitting a play button in GUI.
> By typing in plot.data in the command, I am seeing the numbers only.
> It seems like it is the iterperiod because the number is increased by 100.
> Am I missing something?
> Thanks a lot Jan!
>
> Seungcheol
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Hello Jan,

I just did it and it appeared as "nan" instead of 0 or numbers.
is this meaning my simulation is not working? or something went worng?
I am sorry for bothering you but you helped me a lot.
I sincerely appreciate it.

Seungcheol

Sorry if my suggestions are wrong, I can not check out your script right now.
But I think it is just a problem in order of commands.
If you plot before calculation, you will never get data in your plot.
So this part...

plot.plots={'i':('b0displ')}
plot.plot()
plot.live=True
plot.autozoom=True
O.dt=0.5*PWaveTimeStep()
qt.View()

...should be at the end of your script (after O.run command)

Hello Seungcheol,

yes, it means that something went wrong :-) Try to run the simulation step
by step (so instead of play button you use step button) and see if there is
some warning or errors..

at the beginning of the simulation and after each step (or a few steps) ou
can also check
for b in O.bodies: print b.state.pos
to see where the problems first occur (normally it should be normal
numbers, so when you see NaN, it is a signal of a problem).

Jan

PS: I just checked adding of Vector3 to plot.addData, so meybe the key
would not be 'b0displ' but rather something like 'b0displ_x'..

2013/7/24 Seungcheol Yeom <email address hidden>

> Question #232941 on Yade changed:
> https://answers.launchpad.net/yade/+question/232941
>
> Seungcheol Yeom posted a new comment:
> Hello Jan,
>
> I just did it and it appeared as "nan" instead of 0 or numbers.
> is this meaning my simulation is not working? or something went worng?
> I am sorry for bothering you but you helped me a lot.
> I sincerely appreciate it.
>
> Seungcheol
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Hello Christian,
from this point of view, the Seungcheol's script is ok. This is standard
way how to do live plots (starting plotting before the simulation is run
and updating every while), see e.g.
yade/py/examples/simple-scene/simple-scene-plot.py (with commenting O.run()
command) or many other example scripts.
Jan

2013/7/24 Christian Jakob <email address hidden>

> Question #232941 on Yade changed:
> https://answers.launchpad.net/yade/+question/232941
>
> Christian Jakob posted a new comment:
> Sorry if my suggestions are wrong, I can not check out your script right
> now.
> But I think it is just a problem in order of commands.
> If you plot before calculation, you will never get data in your plot.
> So this part...
>
> plot.plots={'i':('b0displ')}
> plot.plot()
> plot.live=True
> plot.autozoom=True
> O.dt=0.5*PWaveTimeStep()
> qt.View()
>
> ...should be at the end of your script (after O.run command)
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Sorry for my totally dump suggestion in the last post.

I checked out your script today and recognized, that you tried to plot a vector b0displ over iteration number.
Small workaround to fix that:

def addPlotData():
 d0 = O.bodies[0].state.displ()
 d1 = O.bodies[1].state.displ()
 plot.addData(i=O.iter,b0displ=d0[2],b1displ=d1[2])

But you will get zero all the time, because meniscii is not created. So this you can do by changing these lines:

r = 1e-4 #particle radius
h = 1e-6 #praticle distance

#create two sphere paticles#
O.bodies.append([
   utils.sphere(center=(0,0,0),radius=r,fixed=False),
   utils.sphere((0,0,2*r-h),r)
])

hih,

c

Thanks a lot Jan and Christian for your help.
It is working and now I am changing the configuration to check the capillary pressure.
I sincerely appreciate it.
I might ask it again if I have any questions.

Seungcheol

Hello Jan,

I am running the simulation but it seems that it does not influence the particle movement.
When I read (Law2_ScGeom_CapillaryPhys_Capillarity), it sounds like the meniscii is automatically generated based on the capillary pressure and the particle configuration.
So, I tried to make a high capillary pressure which is leading to a small radius of meniscus as well as a low capillary pressure leading to big radius of meniscus.
None of them don't affect the movement of particles.
In addition, I have also changed the particle distant and it did not work.
Any suggestions?
I appreciate your help.

Seungcheol

I am sorry I gorfot to write to thank you for Christian and I need your help.

... , where
aabb = axis aligned bounding box:

https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.Aabb

Thanks Christian Jakob, that solved my question.

Hello Christian,

OMG, you made my simulation work. Thank you very much!
By the way, I added aabbEnlargeFactor greater than 1 as well as interactionDetectionFactor greater than 1.
Can you please let me know what that mean is if you are available?

Seungcheol

A problem is that liquid bridges are stored inside contact data, so normaly there is no bridge if no contact exist.
A solution is to force yade to generate interactions even when there is no contact. As soon as distance is less (radius1+radius2)* interactionDetectionFactor an interaction is created and updated.
You have to be consistent in the interaction detection and contact update, this is why you have the two factors in two different classes.

As far as I know one of condition for liquid bridge is to have a
"normal" (e.g. penetration) contact before capillar force starts to act.
In this case you need only not to delete the interaction after the
distance between particles increases.

Anton

@Anton
That is right, but not deleting is not enough. You need to also create them at one point (hence the increased aabb).
Moreover, you need to also update the normal, which would not be done if we were just "not deleting".