an issue with HertzianViscoElasticFrictionPrms interaction

Hi Horacio,

There is no problem to use "RotSphere" particles with the HertzianViscoElasticFrictionPrms interaction group. The problem you are experiencing arises from the use of the Checkpoint-Restart facility. Checkpoint-Restart is designed only for recommencing simulations that use identical interaction groups. It is not possible to change the interaction groups in the restart script from say RotFrictionPrms to HertzianViscoElasticFrictionPrms. The reason is that the restart file contains information about each of the RotFriction pair-interactions so as to set up those interactions in exactly the same state as before the restart. For RotFrictionPrms, the total shear displacement and current slipping state must be stored in the Checkpoint-Restart file.

One possible way to avoid this issue is to use HertzianViscoElasticFrictionPrms when preparing your system of grains, prior to saving the state to a restart file. As long as the original script to prepare your system of grains has identical interaction groups as the restart script, there shouldn't be a problem.

If you really need to use RotFrictionPrms when preparing the system of grains, you could convert the data in the restart file into an input geometry file. The format of the geometry file is explained here:
https://answers.launchpad.net/esys-particle/+faq/877

The downside of using an input geometry file is that you lose information on the previous positions, velocities and orientations of the particles. Without knowing more about your particular application I cannot say whether this is likely to be a problem or not.

I hope these suggestions help and have fun!

Dion.

Dear Dion,
Thanks for your reply and sorry for not answering sooner to thank you. Everything you say on your answer makes sense and I've tried it but I still get problems using these Interactions when reading from a Checkpoint-Restart file. I haven't had time to fully test implementing the initial setup using a geometry file but I will in the next days and hopefully I'll be able to use the HertzianViscoElasticFrictionPrms for the simulation, since for the preparation is not that important what interactions I use. Just for completeness, I have tried different "combinations" of interactions and I still get errors, for instance using LinearDashpotPrms as in the following instructions:

sim.createInteractionGroup(
   NRotElasticPrms(
      name = "elasticSpheres",
      normalK = elasticContactModulus
   )
)
sim.createInteractionGroup(
   LinearDashpotPrms(
      name = "my_dashpots",
      damp = frictionCoefficient,
      cutoff = 1.2
   )
)
I can not use the Checkpoint-Restart file, I get the same message as before: "number of dynamic interaction groups differ between snapshot and script!" with no success. As I say I will work on using a geometry file, probably the best way to set this problem, and I will let you know if I either I'm successful or if I need extra directions. In the meantime thank you very much for your help!
-Horacio

Dear Dion,
Thanks for your reply and sorry for not answering sooner to thank you.
Everything you say on your answer makes sense and I've tried it but I still
get problems using these Interactions when reading from a Checkpoint-Restart
file.

I am now using a dump checkpoint file into an input geometry file for the
problem I am working with (simple collapse of a granular block) and the
interactions seem to work alright. I will let you know about any other
issues in the near future. Meantime thanks for your help!
Horacio

On Fri, Jun 25, 2010 at 3:01 PM, Dion Weatherley <
<email address hidden>> wrote:

> Your question #115759 on ESyS-Particle changed:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> Status: Open => Answered
>
> Dion Weatherley proposed the following answer:
> Hi Horacio,
>
> There is no problem to use "RotSphere" particles with the
> HertzianViscoElasticFrictionPrms interaction group. The problem you are
> experiencing arises from the use of the Checkpoint-Restart facility.
> Checkpoint-Restart is designed only for recommencing simulations that
> use identical interaction groups. It is not possible to change the
> interaction groups in the restart script from say RotFrictionPrms to
> HertzianViscoElasticFrictionPrms. The reason is that the restart file
> contains information about each of the RotFriction pair-interactions so
> as to set up those interactions in exactly the same state as before the
> restart. For RotFrictionPrms, the total shear displacement and current
> slipping state must be stored in the Checkpoint-Restart file.
>
> One possible way to avoid this issue is to use
> HertzianViscoElasticFrictionPrms when preparing your system of grains,
> prior to saving the state to a restart file. As long as the original
> script to prepare your system of grains has identical interaction groups
> as the restart script, there shouldn't be a problem.
>
> If you really need to use RotFrictionPrms when preparing the system of
> grains, you could convert the data in the restart file into an input
> geometry file. The format of the geometry file is explained here:
> https://answers.launchpad.net/esys-particle/+faq/877
>
> The downside of using an input geometry file is that you lose
> information on the previous positions, velocities and orientations of
> the particles. Without knowing more about your particular application I
> cannot say whether this is likely to be a problem or not.
>
> I hope these suggestions help and have fun!
>
> Dion.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
>
> https://answers.launchpad.net/esys-particle/+question/115759/+confirm?answer_id=0
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> You received this question notification because you are a direct
> subscriber of the question.
>

--
M. Sc. Horacio Tapia McClung
Estudiante de Doctorado
Posgrado de Ingenieria Mecanica
Facultad de Ingenieria/Laboratorio de Fluidos Complejos, IIM
UNAM
+1(831) 824-4039
--
Our virtues and our failings are inseparable, like force and matter. When
they separate, man is no more. -Nikola Tesla

Hi Horacio,

Great to hear you managed to find a solution. To solve your checkpoint-restart problem, I would need to look closely at your simulation scripts. Checkpoint-restart is still considered "beta" by the developers but some users have used it quite successfully. Perhaps there are some bugs in checkpoint-restart for less common combinations of interaction groups?

Have fun with ESyS-Particle and I look forward to hearing how you progress.

Cheers,

Dion.

Hello again Dion
Thanks for your help. Here is an update of the situation I am working with. I managed to solve the initial setup considerations as I described on my previous message and tried to run the simulation using a HertzianViscoElasticFriction interaction group as below, but the simulation does not happen, after reading the .geo file the program only "creates" at most 3 particles and ends. Whenever I change the interaction group for a LinearDashpot or simply remove the Hertzian interaction, the simulation runs fine. Do you know any issue with this interaction group that needs to be taken into account? The same happens when I use a HertzianViscoElastic interaction so it makes me think there is something with the Hertzian interactions.
Thanks for your help and attention,
Horacio

sim.createInteractionGroup (
   HertzianViscoElasticFrictionPrms (
      name = "HertzViscousSpheres",
      A=0.01,
      E=1.0e9,
      nu=0.3,
      dynamicMu = 0.4, #frictionCoefficient,
      shearK = elasticContactModulus
   )
)

Hello again Dion

Thanks for your help. Here is an update of the situation I am working with.
I managed to solve the initial setup considerations as I described on my
previous message and tried to run the simulation using a
HertzianViscoElasticFriction interaction group as below, but the simulation
does not happen, after reading the .geo file the program only "creates" at
most 3 particles and ends. Whenever I change the interaction group for a
LinearDashpot or simply remove the Hertzian interaction, the simulation runs
fine. Do you know any issue with this interaction group that needs to be
taken into account? The same happens when I use a HertzianViscoElastic
interaction so it makes me think there is something with the Hertzian
interactions.
Thanks for your help and attention,
Horacio

sim.createInteractionGroup (
  HertzianViscoElasticFrictionPrms (
     name = "HertzViscousSpheres",
     A=0.01,
     E=1.0e9,
     nu=0.3,
     dynamicMu = 0.4, #frictionCoefficient,
     shearK = elasticContactModulus
  )
)

On Thu, Jul 15, 2010 at 5:17 PM, Dion Weatherley <
<email address hidden>> wrote:

> Your question #115759 on ESyS-Particle changed:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> Status: Open => Answered
>
> Dion Weatherley proposed the following answer:
> Hi Horacio,
>
> Great to hear you managed to find a solution. To solve your checkpoint-
> restart problem, I would need to look closely at your simulation
> scripts. Checkpoint-restart is still considered "beta" by the developers
> but some users have used it quite successfully. Perhaps there are some
> bugs in checkpoint-restart for less common combinations of interaction
> groups?
>
> Have fun with ESyS-Particle and I look forward to hearing how you
> progress.
>
> Cheers,
>
> Dion.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
>
> https://answers.launchpad.net/esys-particle/+question/115759/+confirm?answer_id=3
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> You received this question notification because you are a direct
> subscriber of the question.
>

--
M. Sc. Horacio Tapia McClung
Estudiante de Doctorado
Posgrado de Ingenieria Mecanica
Facultad de Ingenieria/Laboratorio de Fluidos Complejos, IIM
UNAM
+1(831) 824-4039
--
Our virtues and our failings are inseparable, like force and matter. When
they separate, man is no more. -Nikola Tesla

Hi Horatio,

I'm not sure what is the cause of your problems. When you include the Hertzian interaction group, do you remove the NRotElastic and the LinearDashpot interaction groups or only the LinearDashpot interaction group? I would remove both because the Hertzian IG includes elastic repulsion.

Perhaps you could post your complete scripts in your reply or if you prefer, email them to me and I'll take a look.

Cheers,

Dion.

Hi Dion
Thanks again for your reply. I usually run the simulation with either one of
LinearDashpot or NRotElastic plus Hertzian IG but I am trying your
suggestion as i write this message and the result is the same. Any time I
use the Hertzian IG the code does not work well. As you suggest, I am
posting the whole script (excluding the .geo and parameters files) for you
to take a look and if you notice anything that could point me in the right
direction I will greatly appreciate. Please tell me if you would like the
complete setup so that you can see what I am doing at each step. In the
meantime here is the part that where I want to use the Hertzian IG
Many thanks again, and best regards
Horacio

%%%% Code follows below this line

# import the appropriate ESyS-Particle modules:
from esys.lsm import *
from esys.lsm.util import *
from esys.lsm.geometry import *
from columnParams import *
from math import sqrt
import os

# read parameters from file. This file is created separately as described on
the documentation
filename="column.prm"
params = Parameters(filename)
# make local copies of relevant parameters:
#minGrainRadius = params.minRadius
#maxGrainRadius = params.maxRadius
elasticContactModulus = params.elasticContactModulus
contactStrength = params.contactStrength
frictionCoefficient = params.frictionCoefficient
viscosity = params.viscosity
poissonRatio = params.poissonRatio
columnWidth = params.columnWidth
columnHeight = params.columnHeight
dt = params.dt
numTimeSteps = params.numTimeSteps
numSubdomainsX = params.numSubdomainsX
numSubdomainsY = params.numSubdomainsY

#instantiate a simulation object and
#initialise the neighbour search algorithm:
sim = LsmMpi (
      numWorkerProcesses = numSubdomainsX * numSubdomainsY,
      mpiDimList = [numSubdomainsX, numSubdomainsY, 1])
#sim.initNeighbourSearch (
sim.initVerletModel (
   particleType = "RotSphere",
   gridSpacing = 2.5000,
   verletDist = 0.1000
)
sim.force2dComputations(True)

nnumTimeSteps = 6*numTimeSteps
#specify the number of timesteps and the timestep increment:
sim.setNumTimeSteps (nnumTimeSteps)
sim.setTimeStepSize (dt)

# read geometry file
column_geo = "column.geo"
#print column_geo
sim.readGeometry(column_geo)

#add a wall as a floor for the model:
sim.createWall (
   name = "floor",
   posn = Vec3(0.0,0.0,0.0),
   normal = Vec3(0.0,1.0,0.0)
)
#add a left side wall to the model:
sim.createWall (
   name = "left_wall",
   posn = Vec3(0.0,0.0,0.0),
   normal = Vec3(1.0,0.0,0.0)
)
%%%%%%%%% below are the interactions I have set.
#specify that particles undergo frictional interactions:
#sim.createInteractionGroup (
# RotFrictionPrms (
# name = "friction",
# normalK = elasticContactModulus,
# dynamicMu = 0.4, #frictionCoefficient,
# staticMu = frictionCoefficient,
# shearK = elasticContactModulus,
# scaling = True
# )
#)
#specify the type of interactions between colliding particles:
#sim.createInteractionGroup(
# NRotElasticPrms(
# name = "elasticSpheres",
# normalK = elasticContactModulus
# )
#)
#sim.createInteractionGroup(
# LinearDashpotPrms(
# name = "my_dashpots",
# damp = frictionCoefficient,
# cutoff = 1.2
# )
#)

# specify a Hertzian contact interaction groupup
sim.createInteractionGroup (
   HertzianViscoElasticFrictionPrms (
      name = "HertzViscousSpheres",
      A=0.01,
      E=1.0e9,
      nu=0.3,
      dynamicMu = 0.4, #frictionCoefficient,
      shearK = elasticContactModulus
   )
)
# Interactions with walls
#specify that particles undergo elastic repulsion
#from the floor :
sim.createInteractionGroup (
   NRotElasticWallPrms (
      name = "fl_repel",
      wallName = "floor",
      normalK = elasticContactModulus
   )
)
#specify that particles undergo elastic repulsion
#from the left side wall:
sim.createInteractionGroup (
   NRotElasticWallPrms (
      name = "lw_repel",
      wallName = "left_wall",
      normalK = elasticContactModulus
   )
)
# specify elastic repulsion of particles with right wall
#r_repel = NRotElasticWallPrms("rw_repel",
"right_wall",0.0*elasticContactModulus)
#sim.createInteractionGroup (r_repel)

#specify the direction and magnitude of gravity:
sim.createInteractionGroup (
   GravityPrms (
      name = "gravity",
      acceleration = Vec3(0.0000, -9.8100, 0.0000)
   )
)
#add bulk viscosity to damp particle oscillations:
sim.createInteractionGroup (
   LinDampingPrms (
      name = "lin_bulk_viscosity",
      viscosity = viscosity,
      maxIterations = 50
   )
)
sim.createInteractionGroup (
   RotDampingPrms (
      name = "rot_bulk_viscosity",
      viscosity = viscosity,
      maxIterations = 50
   )
)

#add a CheckPointer to store simulation data:
sim.createCheckPointer (
   CheckPointPrms (
      fileNamePrefix = fileNamePrefix,
      beginTimeStep = numTimeSteps,
      endTimeStep = nnumTimeSteps,
      timeStepIncr = 1000
   )
)

#execute the simulation:
sim.run()

On Sat, Jul 31, 2010 at 5:18 PM, Dion Weatherley <
<email address hidden>> wrote:

> Your question #115759 on ESyS-Particle changed:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> Status: Open => Answered
>
> Dion Weatherley proposed the following answer:
> Hi Horatio,
>
> I'm not sure what is the cause of your problems. When you include the
> Hertzian interaction group, do you remove the NRotElastic and the
> LinearDashpot interaction groups or only the LinearDashpot interaction
> group? I would remove both because the Hertzian IG includes elastic
> repulsion.
>
> Perhaps you could post your complete scripts in your reply or if you
> prefer, email them to me and I'll take a look.
>
> Cheers,
>
> Dion.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
>
> https://answers.launchpad.net/esys-particle/+question/115759/+confirm?answer_id=6
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> You received this question notification because you are a direct
> subscriber of the question.
>

--
M. Sc. Horacio Tapia McClung
Estudiante de Doctorado
Posgrado de Ingenieria Mecanica
Facultad de Ingenieria/Laboratorio de Fluidos Complejos, IIM
UNAM
+1(831) 824-4039
--
Our virtues and our failings are inseparable, like force and matter. When
they separate, man is no more. -Nikola Tesla

Hi Horacio,

Could you please also post a list of the parameters you are using (in your "column.prm" file)? It would also help if you could describe the initial geometry you are using but I doubt that is the source of the problem.

Cheers,

Dion.

Hello Dion
I am sorry, of course I can post the information you request, I am sorry for
not doing it. Here I post all the information you need to try and run a
simulation with the code I sent you before. As always thanks for your help!
Horacio

########## Sample parameters file (note that some of the parameters listed
on this file are not used on the program)

column.prm

elasticContactModulus 100000.000
contactStrength 100.000
frictionCoefficient 0.100
viscosity 0.100
poissonRatio 0.300
dt 1.000e-04
columnWidth 19.000
columnHeight 7.000
plateWidth 0.000
minRadius 0.500
maxRadius 0.600
minGrainRadius 2.000
maxGrainRadius 10.000
numTimeSteps 100000
numSubdomainsX 1
numSubdomainsY 1

########## Geometry file of a simple 2D step or column. bBelow is a sample
with 93 particles.

column.geo:

LSMGeometry 1.2
BoundingBox 0 0 0.0 20 7.5 0.0
PeriodicBoundaries 0 0 0
Dimension 2D
BeginParticles
Simple
93
0.567232 0.56728 0 0.567344 215 123654
1.67498 0.540893 0 0.540936 0 123654
0.580736 1.71528 0 0.580805 514 123654
1.54508 2.37556 0 0.58801 1 123654
2.14769 1.46903 0 0.500697 4 123654
1.36435 3.9075 0 0.521387 2 123654
0.598404 3.0906 0 0.598445 127 123654
2.16854 3.27173 0 0.503765 6 123654
2.81018 0.595579 0 0.595736 110 123654
3.95038 0.545632 0 0.545684 149 123654
3.42516 1.5077 0 0.504553 3 123654
3.87522 2.46083 0 0.549627 9 123654
2.75734 2.34383 0 0.56557 5 123654
4.46383 1.52542 0 0.534347 8 123654
4.98689 2.48943 0 0.562487 14 123654
2.50419 4.29775 0 0.575794 7 123654
3.24902 3.4072 0 0.585218 11 123654
4.39327 3.44155 0 0.559563 10 123654
3.81688 4.43089 0 0.585455 12 123654
5.08411 0.588891 0 0.589017 331 123654
6.25714 0.584082 0 0.584184 109 123654
5.60129 1.58356 0 0.53218 13 123654
6.68747 1.6392 0 0.555516 18 123654
7.3506 0.511666 0 0.511719 203 123654
7.34443 2.49391 0 0.522646 24 123654
5.03221 4.32095 0 0.527526 16 123654
6.07445 2.52997 0 0.525905 19 123654
5.5325 3.40241 0 0.501163 15 123654
6.73701 3.3836 0 0.554695 20 123654
6.10101 4.2772 0 0.542177 21 123654
7.2787 4.39809 0 0.595425 26 123654
5.35728 5.33479 0 0.537192 17 123654
6.43829 5.3097 0 0.544077 22 123654
7.8111 1.4982 0 0.577085 23 123654
9.53295 0.578416 0 0.578599 74 123654
8.40831 0.546511 0 0.546582 105 123654
8.95212 1.56457 0 0.565978 28 123654
8.37874 2.47719 0 0.511854 29 123654
9.5268 2.47595 0 0.511469 34 123654
7.88449 3.42326 0 0.552319 25 123654
8.95414 3.33494 0 0.520976 30 123654
8.38117 4.35877 0 0.506915 31 123654
9.39178 4.2661 0 0.507935 36 123654
7.55271 5.5597 0 0.598067 27 123654
8.74063 5.40554 0 0.599866 32 123654
9.83906 5.21362 0 0.515235 37 123654
10.6324 0.520293 0 0.520305 450 123654
10.0604 1.56727 0 0.542364 33 123654
11.7219 0.560378 0 0.560524 188 123654
11.1789 1.55893 0 0.576359 38 123654
12.3054 1.57447 0 0.550319 43 123654
10.5931 2.5119 0 0.542283 39 123654
10.0007 3.3702 0 0.500564 35 123654
10.4325 4.34607 0 0.535923 41 123654
11.7547 2.54786 0 0.568104 44 123654
11.1193 3.48586 0 0.564863 40 123654
11.5869 4.4995 0 0.551484 46 123654
12.2698 3.58496 0 0.58999 45 123654
10.8771 5.31216 0 0.527501 42 123654
11.9679 5.58147 0 0.595631 47 123654
12.869 0.58683 0 0.586943 77 123654
13.4298 1.54937 0 0.527183 48 123654
12.8672 2.47019 0 0.507001 49 123654
14.0118 0.556288 0 0.556367 79 123654
14.5226 1.5551 0 0.565615 53 123654
13.371 3.43548 0 0.521537 50 123654
12.713 4.6697 0 0.581879 51 123654
13.9323 2.51222 0 0.558966 54 123654
14.4591 3.54226 0 0.571826 55 123654
13.7772 4.39229 0 0.517936 56 123654
14.9976 4.57263 0 0.590871 61 123654
13.1368 5.71832 0 0.549179 52 123654
14.1881 5.38274 0 0.554391 57 123654
15.1037 0.535745 0 0.535764 171 123654
16.2099 0.570782 0 0.571035 147 123654
15.658 1.56898 0 0.569894 58 123654
16.7983 1.5437 0 0.566113 63 123654
17.3332 0.552438 0 0.552513 347 123654
15.0739 2.56378 0 0.583915 59 123654
16.229 2.51271 0 0.533235 64 123654
15.7432 3.49901 0 0.566246 60 123654
16.1549 4.5544 0 0.566678 66 123654
16.8107 3.45111 0 0.502301 65 123654
17.3051 2.53012 0 0.543038 69 123654
17.3 4.4599 0 0.582402 71 123654
15.6167 5.57913 0 0.590831 62 123654
16.7909 5.46148 0 0.541171 67 123654
17.8795 1.47236 0 0.517537 68 123654
18.4429 0.557144 0 0.557205 377 123654
18.4216 2.42478 0 0.57848 197 123654
17.9211 3.47326 0 0.583455 70 123654
18.441 4.4905 0 0.559036 141 123654
17.9104 5.49716 0 0.578962 72 123654
EndParticles
BeginConnect
0
EndConnect

I believe that is all the input information you need to try and reproduce
the problem I get. If you have any other request I will be glad to send you
more information as I am really looking forward to continue using esys as a
research tool :)

Thanks again and again :)
Horacio

On Thu, Aug 5, 2010 at 6:05 PM, Dion Weatherley <
<email address hidden>> wrote:

> Your question #115759 on ESyS-Particle changed:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> Status: Open => Needs information
>
> Dion Weatherley requested for more information:
> Hi Horacio,
>
> Could you please also post a list of the parameters you are using (in
> your "column.prm" file)? It would also help if you could describe the
> initial geometry you are using but I doubt that is the source of the
> problem.
>
> Cheers,
>
> Dion.
>
> --
> To answer this request for more information, you can either reply to
> this email or enter your reply at the following page:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> You received this question notification because you are a direct
> subscriber of the question.
>

--
M. Sc. Horacio Tapia McClung
Estudiante de Doctorado
Posgrado de Ingenieria Mecanica
Facultad de Ingenieria/Laboratorio de Fluidos Complejos, IIM
UNAM
+1(831) 824-4039
--
Our virtues and our failings are inseparable, like force and matter. When
they separate, man is no more. -Nikola Tesla

This question was expired because it remained in the 'Open' state without activity for the last 15 days.

Question not yet answered

Hi Dion
Thanks for reopening this question, I still don't know what could be the
issue. I actually wonder if you received the information I sent you (the
code and initial geometry file).
I will be checking and continue experimenting
Thanks again
Horacio

On Wed, Aug 25, 2010 at 4:46 PM, Dion Weatherley <
<email address hidden>> wrote:

> Your question #115759 on ESyS-Particle changed:
> https://answers.launchpad.net/esys-particle/+question/115759
>
> Status: Expired => Open
>
> Dion Weatherley changed the question status:
> Question not yet answered
>
> --
> You received this question notification because you are a direct
> subscriber of the question.
>

--
M. Sc. Horacio Tapia McClung
Estudiante de Doctorado
Posgrado de Ingenieria Mecanica
Facultad de Ingenieria/Laboratorio de Fluidos Complejos, IIM
UNAM
+1(831) 824-4039
--
Our virtues and our failings are inseparable, like force and matter. When
they separate, man is no more. -Nikola Tesla

Hi Horatio,

Yes, I received your code and parameter files. Unfortunately, I have been too busy with my own work for the past two weeks to investigate your problem further. Hopefully I'll get to it in the next few days.

Cheers,

Dion.