numerical instability and small time steps

Theoreticaly, the timesteps should be the same. What contact laws are you using and are you sure the contact parameters are identical?
Did you try setting the timestep update interval to 100 in PFC as well? I don't remember the fish command but I know it can be changed.
Could you send the script?

On 21 October 2011 17:25, Chareyre <email address hidden>wrote:

> Question #175394 on Yade changed:
> https://answers.launchpad.net/yade/+question/175394
>
> Status: Open => Answered
>
> Chareyre proposed the following answer:
> Theoreticaly, the timesteps should be the same. What contact laws are you
> using and are you sure the contact parameters are identical?
>
I remember you were using the HertzMindlin law. If that is still the case,
remember that kn and ks are updated every time step and that is accounted
for in GlobalStiffnessTimeStepper.

Chiara

> Did you try setting the timestep update interval to 100 in PFC as well? I
> don't remember the fish command but I know it can be changed.
> Could you send the script?
>
> --
> 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
>

I am using hertz model with following properties:

poisson_ratio = 0.3
young_modulus = 2*shear_modulus*(1+poisson_ratio)
friction_coeff = 2
angle = math.atan(friction_coeff)
rho_p = 2650

My materials and engines:

id_FacetMat=O.materials.append(FrictMat(young=young_modulus,poisson=poisson_ratio,frictionAngle=angle))
id_SphereMat=O.materials.append(FrictMat(young=young_modulus,poisson=poisson_ratio,density=rho_p,frictionAngle=angle))

FacetMat=O.materials[id_FacetMat]
SphereMat=O.materials[id_SphereMat]

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Facet_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_CapillaryPhys(label='ContactModel')],
  [Law2_ScGeom_FrictPhys_CundallStrack()],
  ),
 Law2_ScGeom_CapillaryPhys_Capillarity(CapillaryPressure=10000),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=5,timestepSafetyCoefficient=0.8, defaultDt=4*utils.PWaveTimeStep()),
 GravityEngine(gravity=(0,0,-9.81)),
 NewtonIntegrator(damping=local_damping),
]
ContactModel.betan=viscous_normal # normal direction
ContactModel.betas=viscous_shear # shear direction
ContactModel.useDamping=True

I think this is ok, or did I do a mistake?

In PFC I use the same properties shear mod.= 3e7, poisson = 0.3 and rho = 2650.
The only difference is the wall, where I use linear model with stiffnesses kn=1e10 and ks=0 in PFC.

Can I do the same in Yade? ... Using hertz model for particles and linear model for the facets?

christian

p.s. with timeStepUpdateInterval=5 calculation is stable ... PFC uses 1 as default update interval, I did not try to set it to 100 ...

If you use Ip2_FrictMat_FrictMat_CapillaryPhys + Law2_ScGeom_FrictPhys_CundallStrack, then you are not simulating Hertz-Mindlin, but a linear law instead.

These lines have no effect:
ContactModel.betan=viscous_normal # normal direction
ContactModel.betas=viscous_shear # shear direction
ContactModel.useDamping=True

The "contact model" (physics functor) has no such attributes (see https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.Ip2_FrictMat_FrictMat_CapillaryPhys). It has in fact no attributes at all, so you are just adding data to it (python allows that) but they are never used in computations.

I think it explains the differences in time-steps since contact stiffness are most likely very different (the straight way to check that is to compare directly contacts kn/ks in Yade and PFC - for Yade use simulation inspector to see the values).

In order to simulate Hertz contact + capillary forces, I think it will need an additional Ip2 functor implemented, that will generate the capillary data in addition to what is already defined by Ip2_FrictMat_FrictMat_MindlinPhys. There is nothing to implement in terms of physics, it's only a matter of declaring things correctly and creating the good Ip type.

For different contact laws with walls, it is theoreticaly possible, but I think that, again, it will need a little bit of work on the c++ side.

Currently, it would work if combining functors and laws like this for instance:
InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_Dem3Dof()],
  [Ip2_FrictMat_FrictMat_CapillaryPhys(label='ContactModel')],
  [Law2_ScGeom_FrictPhys_CundallStrack(), Law2_Dem3Dof_FrictPhys_CundallStrack()],

As you see here, sphere-sphere and sphere-facet will have different geometry types (ScGeom vs. Dem3Dof), hence they will be dispatched to different contact laws (even if the two contact laws I mention here are doing the same thing actually).
In your case, it would not work unfortunately, because you need ScGeom for both contact groups. A solution would be to make one contact group artificially different by deriving another geometry type from ScGeom, that would be exactly the same except the name. It is easy to do with inheritance but sounds like a hack.
This is actually an interesting problem since it shows a design problem in Yade.

Last note: you don't need Ig2_Sphere_Sphere_ScGeom6D. The "6D" version is for contact moments. You can use the simplest ScGeom. It is harmless to use 6D though, except for memory usage.

Thank you for your fast and extensive answer.

So the problem is, that I set a Hertz model, but I get a linear model.
But I can not set a linear model with the CapillaryPhys functor...
This is very confusing ...

As you said I have to create a new Ip2 functor to use it for my special problem.
Can you estimate how much time this will cost for a c++ noob? :D

christian

cost =
time to read a c++ book, especially inheritance
+
time to understand yade's dispatching mechanism (if not clear yet)
+
30min typing and compiling if everything goes well.

I may be able to help a little but not before next week.

On 24 October 2011 12:55, Christian Jakob <
<email address hidden>> wrote:

> Question #175394 on Yade changed:
> https://answers.launchpad.net/yade/+question/175394
>
> Christian Jakob posted a new comment:
> Thank you for your fast and extensive answer.
>
> So the problem is, that I set a Hertz model, but I get a linear model.
>
No, you do not set the Hertz model at all (the corresponding law functor
would be:
https://yade-dem.org/doc/yade.wrapper.html?highlight=mindlin#yade.wrapper.Law2_ScGeom_MindlinPhys_Mindlin).
You set a linear model and at the same time you introduce some attributes in
python, like ContactModel.betan, which have no effect. As Bruno said, you
need a new Ip functor if you want to use the existing capillary model with
the Hertzian law. Chiara

> But I can not set a linear model with the CapillaryPhys functor...
> This is very confusing ...
>
> As you said I have to create a new Ip2 functor to use it for my special
> problem.
> Can you estimate how much time this will cost for a c++ noob? :D
>
> christian
>
> --
> 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
>

How should we call the new Ip2 functor?

I suggest: Ip2_FrictMat_FrictMat_MindlinCapillaryPhys

Better ideas?

christian

> I suggest: Ip2_FrictMat_FrictMat_MindlinCapillaryPhys
>
It makes sense.

B

Now, while using Ip2_FrictMat_FrictMat_MindlinCapillaryPhys instead of Ip2_FrictMat_FrictMat_CapillaryPhys, I still get very small time steps (2e-9 s). The problem remains the same ...

This is still an open question for me too. The timesteps should be the same if the contact parameters are really the same (this still needs to be double-checked).
If identical contact parameters and dt give stable integration in PFC but unstable in Yade, it most probably reveals a bug.

I agree. The time step should be the same if computed the SAME way! If Christian is using viscous damping (are you?) in PFC, this will affect the time step determination as well (according to the PFC manual, at least). So the first thing to do is to apply the same time step determination and only then compare the two codes...

It is not double-checked, but quad-checked (maybe even more).

YADE input file:

#GENERAL SETTINGS:
initial_porosity = 0.4
shear_modulus = 3e7
poisson_ratio = 0.3
young_modulus = 2*shear_modulus*(1+poisson_ratio)
friction_coeff = 2
angle = math.atan(friction_coeff)
rho_p = 2650 #density of particles

#SET DAMPING CONSTANTS:
local_damping = 0.01
viscous_normal = 0.021
viscous_shear = 0.8*viscous_normal

PFC input file:

 ;GENERAL SETTINGS:
 initial_porosity = 0.4
 shear_modulus = 3e7
 poisson_ratio = 0.3
 friction_coeff = 2
 rho_p = 2650 ;density of particles

 ;SET DAMPING CONSTANTS:
 local_damping = 0.01
 viscous_normal = 0.021
 viscous_shear = 0.8*viscous_normal

The only thing I could imagine, that is wrong, is

young_modulus = 2*shear_modulus*(1+poisson_ratio)

but wikipedia tells me this is ok...

time step determination in PFC can be different, than in Yade, i dont know.
but i know the time step width, that is stable in PFC (2e-7 s). i can set the same time step in yade, but it is unstable.
GSTS in yade tells me, that the time step should be about 2e-9 s for a stable calculation.

@bruno:
stable = !unstable

... very good^^

On 29 November 2011 09:25, Christian Jakob <
<email address hidden>> wrote:

> Question #175394 on Yade changed:
> https://answers.launchpad.net/yade/+question/175394
>
> Christian Jakob posted a new comment:
> It is not double-checked, but quad-checked (maybe even more).
>
Did you check the computation of the time step in PFC in case of viscous
damping and Hertzian law? This has to be quad-checked. In Yade there is no
such a computation of time step that considers viscous damping at the
contact (AFAIK).

C.

>
> YADE input file:
>
> #GENERAL SETTINGS:
> initial_porosity = 0.4
> shear_modulus = 3e7
> poisson_ratio = 0.3
> young_modulus = 2*shear_modulus*(1+poisson_ratio)
> friction_coeff = 2
> angle = math.atan(friction_coeff)
> rho_p = 2650 #density of particles
>
> #SET DAMPING CONSTANTS:
> local_damping = 0.01
> viscous_normal = 0.021
> viscous_shear = 0.8*viscous_normal
>
> PFC input file:
>
> ;GENERAL SETTINGS:
> initial_porosity = 0.4
> shear_modulus = 3e7
> poisson_ratio = 0.3
> friction_coeff = 2
> rho_p = 2650 ;density of particles
>
> ;SET DAMPING CONSTANTS:
> local_damping = 0.01
> viscous_normal = 0.021
> viscous_shear = 0.8*viscous_normal
>
> The only thing I could imagine, that is wrong, is
>
> young_modulus = 2*shear_modulus*(1+poisson_ratio)
>
> but wikipedia tells me this is ok...
>
> --
> 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
>

The timestep determination doesn not really matter here. It should be fixed manually to identical values in both codes, then if only one explodes there is a problem.

Christian, in order to double check contacts parameters, looking at input values is not enough. The best way is to simulate a compression on two spheres and see the displacement/force response. That way, you will know if the contact behavior is really the same.

Ans also the velocity/viscous force response...

Ok, i will create a smaller example to figure out the differences.

@chiara

i looked into PFC manual, but i could not find anything reffered to time step determination in the case when viscous damping is active. there is just a definition of critical time step determination with t_{crit} = \sqrt{\frac{m}{k}}

On 29 November 2011 10:05, Christian Jakob <
<email address hidden>> wrote:

> Question #175394 on Yade changed:
> https://answers.launchpad.net/yade/+question/175394
>
> Christian Jakob posted a new comment:
> Ok, i will create a smaller example to figure out the differences.
>
I do remember that you already checked the Hertz law in Yade and PFC (and
the solution was the same). You could you use that same script just adding
the viscous force. Please, make sure you understand what the contact model
is doing in both codes.

>
> @chiara
>
> i looked into PFC manual, but i could not find anything reffered to time
> step determination in the case when viscous damping is active. there is
> just a definition of critical time step determination with t_{crit} =
> \sqrt{\frac{m}{k}}
>
I am 100% sure there is a modification but I do not have the manual with me
now. As Bruno said, it should not matter provided that the contact model
acts the same way. However, it is good for you to understand/know how/why
the time step should be reduced in case of viscous contact damping.

Chiara

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

Anton mentioned that the script was still changing the radii... it is (was?) a problem for sure.
I may also recall that facets have known bugs... I again saw strange movements of spheres touching the facets recently (looks like at equilibrium, then suddendly bounces back).

@chiara

i searched for "timestep" in PFC manual. i really could not find anything referred to viscous contact damping and timestep determination. it is clear to me, that timesteps should be a bit smaller in that case, but not with factor 100! what i can imagine is, that some code details are not explained in the manual...

i tested a script without viscous damping and the timestep is still small compared to PFC timestep.
i also tested a script without any walls or facets, but it has no effect on timesteps.

as mentioned above, i will try to reactivate the "old" test script with two colliding spheres and come back with the news.

for two colliding spheres i get following timesteps:

PFC with viscous damping: 0.919e-4
PFC without viscous damping: 0.832e-4

YADE with viscous damping: 1.127e-4
YADE without viscous damping: 1.127e-4

conclusion:
- timestep determinations are different
- timestep determ. in PFC takes account for viscous damping

But this does not explain the small timesteps I observed.

Could you check that Yade is table for several bouncing in the 2-sphere case with the dt=~1e-4? This at least would be a good thing.
Spheres also need to have the same mass in Yade and PFC.

Can you please tell us if the force-displacement curves at the contact
match between the two codes? It would be good to know that both with and
without damping.

Thanks, Chiara

On 2 December 2011 09:55, Christian Jakob <
<email address hidden>> wrote:

> Question #175394 on Yade changed:
> https://answers.launchpad.net/yade/+question/175394
>
> Christian Jakob posted a new comment:
> for two colliding spheres i get following timesteps:
>
> PFC with viscous damping: 0.919e-4
> PFC without viscous damping: 0.832e-4
>
> YADE with viscous damping: 1.127e-4
> YADE without viscous damping: 1.127e-4
>
> conclusion:
> - timestep determinations are different
> - timestep determ. in PFC takes account for viscous damping
>
> But this does not explain the small timesteps I observed.
>
> --
> 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
>

> Spheres also need to have the same mass in Yade and PFC.

mass PFC: 4.188790204786e-6
mass YADE: 4.18879020479e-6

force-disp. curves:

PFC without viscous damping:
un ................. vs. Fn
------------ ------------
-1.0000e-004 0.0000e+000
-1.0000e-004 1.2778e-001
-9.3899e-005 1.1626e-001
-8.6133e-005 1.0214e-001
-7.6904e-005 8.6172e-002
-6.6440e-005 6.9198e-002
-5.4985e-005 5.2098e-002
-4.2785e-005 3.5758e-002
-3.0071e-005 2.1071e-002
-1.7056e-005 9.0007e-003

YADE without viscous damping:
un ................. vs. Fn
*first step not available*
-9.4559e-05 0.124283508259
-8.9266e-05 0.117490612129
-8.2429e-05 0.107764352662
-7.4222e-05 0.0956237678404
-6.4845e-05 0.0817047043681
-5.4512e-05 0.0667209759367
-4.3443e-05 0.0514267093956
-3.1849e-05 0.036586864111
-1.9927e-05 0.0229667479151

PFC with viscous damping:
un ................. vs. Fn
------------ ------------
-1.0000e-004 0.0000e+000
-1.0000e-004 1.2778e-001
-9.3899e-005 1.1626e-001
-8.6242e-005 1.0234e-001
-7.7253e-005 8.6760e-002
-6.7174e-005 7.0348e-002
-5.6253e-005 5.3910e-002
-4.4732e-005 3.8227e-002
-3.2834e-005 2.4040e-002
-2.0755e-005 1.2082e-002

YADE with viscous damping:
un ................. vs. Fn
*first step not available*
-9.4594e-05 0.121855832544
-8.9402e-05 0.112865153894
-8.2758e-05 0.101313188504
-7.4857e-05 0.0878212559716
-6.5908e-05 0.0730792438382
-5.6132e-05 0.0578038888246
-4.5744e-05 0.042703264729
-3.4948e-05 0.0284531709081
-2.3928e-05 0.0156941209448

With the method suggested by Anton (replace balls instead of multiply the radius) I get a bigger timestep. But it is still too small: 8e-9 s (compared to PFC: 2e-7 s). So I tried to test, if a fixed timestep of 2e-7 s would be stable:

Yade [2]: *** glibc detected *** /usr/bin/python: malloc(): memory corruption: 0x0000000006e1b7e0 ***
======= Backtrace: =========
/lib/libc.so.6(+0x71ad6)[0x7f714b44dad6]
/lib/libc.so.6(+0x74b6d)[0x7f714b450b6d]
/lib/libc.so.6(__libc_malloc+0x70)[0x7f714b452930]
/usr/lib/libQtGui.so.4(+0x2c8d26)[0x7f7147e59d26]
/usr/lib/libQtGui.so.4(_ZNK12QPainterPath17toSubpathPolygonsERK10QTransform+0x651)[0x7f7147ef6511]
/usr/lib/libQtGui.so.4(_ZNK12QPainterPath13toFillPolygonERK10QTransform+0x23)[0x7f7147ef7113]
/usr/lib/libQtGui.so.4(_ZNK12QPainterPath13toFillPolygonERK7QMatrix+0x32)[0x7f7147ef7482]
/usr/lib/libQtGui.so.4(+0x400c7b)[0x7f7147f91c7b]
/usr/lib/libQtGui.so.4(+0x34bb53)[0x7f7147edcb53]
/usr/lib/libQtGui.so.4(+0x34bc32)[0x7f7147edcc32]
/usr/lib/libQtGui.so.4(_ZN8QPainter13setClipRegionERK7QRegionN2Qt13ClipOperationE+0x369)[0x7f7147ee1dd9]
/usr/lib/kde4/plugins/styles/oxygen.so(+0x5de0e)[0x7f7123fbfe0e]
/usr/lib/libQtCore.so.4(_ZN23QCoreApplicationPrivate29sendThroughObjectEventFiltersEP7QObjectP6QEvent+0x87)[0x7f71478784b7]
/usr/lib/libQtGui.so.4(_ZN19QApplicationPrivate13notify_helperEP7QObjectP6QEvent+0x7c)[0x7f7147d892fc]
/usr/lib/libQtGui.so.4(_ZN12QApplication6notifyEP7QObjectP6QEvent+0x14b)[0x7f7147d8f80b]
/usr/lib/pymodules/python2.6/PyQt4/QtGui.so(+0x4933cf)[0x7f712c60d3cf]
/usr/lib/libQtCore.so.4(_ZN16QCoreApplication14notifyInternalEP7QObjectP6QEvent+0x8c)[0x7f714787909c]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate10drawWidgetEP12QPaintDeviceRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x3bd)[0x7f7147de780d]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate22paintSiblingsRecursiveEP12QPaintDeviceRK5QListIP7QObjectEiRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x518)[0x7f7147de84a8]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate22paintSiblingsRecursiveEP12QPaintDeviceRK5QListIP7QObjectEiRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x349)[0x7f7147de82d9]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate10drawWidgetEP12QPaintDeviceRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x11a)[0x7f7147de756a]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate22paintSiblingsRecursiveEP12QPaintDeviceRK5QListIP7QObjectEiRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x518)[0x7f7147de84a8]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate10drawWidgetEP12QPaintDeviceRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x11a)[0x7f7147de756a]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate22paintSiblingsRecursiveEP12QPaintDeviceRK5QListIP7QObjectEiRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x518)[0x7f7147de84a8]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate10drawWidgetEP12QPaintDeviceRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x11a)[0x7f7147de756a]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate22paintSiblingsRecursiveEP12QPaintDeviceRK5QListIP7QObjectEiRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x518)[0x7f7147de84a8]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate10drawWidgetEP12QPaintDeviceRK7QRegionRK6QPointiP8QPainterP19QWidgetBackingStore+0x11a)[0x7f7147de756a]
/usr/lib/libQtGui.so.4(+0x410a05)[0x7f7147fa1a05]
/usr/lib/libQtGui.so.4(_ZN14QWidgetPrivate16syncBackingStoreEv+0x80)[0x7f7147dd9170]
/usr/lib/libQtGui.so.4(_ZN7QWidget5eventEP6QEvent+0xcd5)[0x7f7147ddf915]
/usr/lib/pymodules/python2.6/PyQt4/QtGui.so(+0x4755e5)[0x7f712c5ef5e5]
/usr/lib/libQtGui.so.4(_ZN19QApplicationPrivate13notify_helperEP7QObjectP6QEvent+0xac)[0x7f7147d8932c]
/usr/lib/libQtGui.so.4(_ZN12QApplication6notifyEP7QObjectP6QEvent+0x14b)[0x7f7147d8f80b]
/usr/lib/pymodules/python2.6/PyQt4/QtGui.so(+0x4933cf)[0x7f712c60d3cf]
/usr/lib/libQtCore.so.4(_ZN16QCoreApplication14notifyInternalEP7QObjectP6QEvent+0x8c)[0x7f714787909c]
/usr/lib/libQtCore.so.4(_ZN23QCoreApplicationPrivate16sendPostedEventsEP7QObjectiP11QThreadData+0x2a4)[0x7f714787c744]
/usr/lib/libQtCore.so.4(+0x1a2b73)[0x7f71478a2b73]
/lib/libglib-2.0.so.0(g_main_context_dispatch+0x1f2)[0x7f71430806f2]
/lib/libglib-2.0.so.0(+0x42568)[0x7f7143084568]
/lib/libglib-2.0.so.0(g_main_context_iteration+0x6c)[0x7f714308471c]
/usr/lib/libQtCore.so.4(_ZN20QEventDispatcherGlib13processEventsE6QFlagsIN10QEventLoop17ProcessEventsFlagEE+0x73)[0x7f71478a26b3]
/usr/lib/libQtGui.so.4(+0x2a819e)[0x7f7147e3919e]
/usr/lib/libQtCore.so.4(_ZN10QEventLoop13processEventsE6QFlagsINS_17ProcessEventsFlagEE+0x32)[0x7f71478779c2]
/usr/lib/libQtCore.so.4(_ZN10QEventLoop4execE6QFlagsINS_17ProcessEventsFlagEE+0xdc)[0x7f7147877d9c]
/usr/lib/libQtCore.so.4(_ZN16QCoreApplication4execEv+0xbb)[0x7f714787ca2b] 09:25:46
/usr/lib/pymodules/python2.6/PyQt4/QtCore.so(+0x98ee5)[0x7f712dce2ee5]
/usr/lib/python2.6/lib-dynload/readline.so(+0x3991)[0x7f7131322991]
/usr/bin/python(PyOS_Readline+0xf2)[0x526a52]
/usr/bin/python[0x49da86]
/usr/bin/python(PyEval_EvalFrameEx+0x5165)[0x4a7ba5]
/usr/bin/python(PyEval_EvalCodeEx+0x911)[0x4a95c1]
/usr/bin/python(PyEval_EvalFrameEx+0x4d12)[0x4a7752]
/usr/bin/python(PyEval_EvalCodeEx+0x911)[0x4a95c1]
/usr/bin/python(PyEval_EvalFrameEx+0x4d12)[0x4a7752]
/usr/bin/python(PyEval_EvalCodeEx+0x911)[0x4a95c1]
/usr/bin/python(PyEval_EvalFrameEx+0x4d12)[0x4a7752]
/usr/bin/python(PyEval_EvalCodeEx+0x911)[0x4a95c1]
/usr/bin/python[0x538a10]
/usr/bin/python(PyObject_Call+0x47)[0x41ef47]
/usr/bin/python[0x427c1f]
/usr/bin/python(PyObject_Call+0x47)[0x41ef47]
/usr/bin/python[0x426edd]
/usr/bin/python(PyObject_Call+0x47)[0x41ef47]
======= Memory map: ========
00400000-0061d000 r-xp 00000000 08:01 2891650 /usr/bin/python2.6
0081d000-0087f000 rw-p 0021d000 08:01 2891650 /usr/bin/python2.6
0087f000-0088e000 rw-p 00000000 00:00 0
0279a000-0858b000 rw-p 00000000 00:00 0 [heap]
7f711b1f6000-7f711b2f6000 rw-s 34363000 00:05 6842 /dev/ati/card0
7f711b2f6000-7f711b3f6000 rw-s 34362000 00:05 6842 /dev/ati/card0
7f711b3f6000-7f711b4f6000 rw-s 34361000 00:05 6842 /dev/ati/card0
7f711b4f6000-7f711b5f6000 rw-s 34360000 00:05 6842 /dev/ati/card0
7f711b5f6000-7f711b6f6000 rw-s 34365000 00:05 6842 /dev/ati/card0
7f711b6f6000-7f711b7f6000 rw-s 3435d000 00:05 6842 /dev/ati/card0
7f711b7f6000-7f711b8f6000 rw-s 3435c000 00:05 6842 /dev/ati/card0
7f711b8f6000-7f711b9f6000 rw-s 34364000 00:05 6842 /dev/ati/card0
7f711baf6000-7f711bbf6000 rw-s 34359000 00:05 6842 /dev/ati/card0
7f711bbf6000-7f711bcf6000 rw-s 34358000 00:05 6842 /dev/ati/card0
Abgebrochen

???

> 7f711b1f6000-7f711b2f6000 rw-s 34363000 00:05 6842 /dev/ati/card0

I think it is a video-driver problem.
Try to start your script without opening a graphic interface.

Anton

From your tables, I conclude that force-displacement realtions are relatively similar (although not exactly equal). It can't explain a factor 100 in dt, provided there is no parameter mismatch introduced in the other scripts (are mass the same in ?).

Another thing is the tangential stiffness, that can also give instabilities. To eliminate this cause, it is enough to set friction=0 and see if it stabilizes the problem.

I don't think the dt problem is connected to viscosity or radii multiplication.
What I'm wondering is if Yade is stable in the simple case of a bouncing sphere (https://www.yade-dem.org/doc/tutorial-examples.html), for the same time-step as PFC. This is a good test: no facets, no change of radii, no viscosity, no shear force, similar force-displacement relations on the normal. If we need different timesteps in this situation, it really narrows down the possible location of the bug.

I dont know what happend, but in my last calculation it seems, that the time step is ok now.
I detected, that I now have 7.25e-7 [s] ... :D
I simply included that isBroken flag into Law2...Capillarity.cpp and everything is fine?!