Interaction of FrictMat and a new class of material

What about the behavior you observed with Yade ? (crash ? with what message if yes ?...)

It would help us to help you I think, as cross reading a bunch of c++ lines is not so easy....

I don't get any error message or crash!

But, the force-displacement law is not following the analytical response and actually the two contacting spheres pass through each other by even a very small force.

When I ask for properties of the interaction like
O.interactions[0,1].phys.kmn

it gives ;
AttributeError: 'NoneType' object has no attribute 'kmn'

but if I ask few steps after, it gives the right answer as:

Yade [4]: O.interactions[0,1].phys.kmn
 -> [4]: 61538461.538461536

>
>
>
> When I ask for properties of the interaction like
> O.interactions[0,1].phys.kmn
>
> it gives ;
> AttributeError: 'NoneType' object has no attribute 'kmn'
>

This means that O.interactions[0,1].phys is None and therefore the
interaction is only potential and not real

Hi,
My first reaction here was "where is the hpp?". It is not shown.
It would be helpful to have the two files. Not just pasted in textbox as indentation and syntax highlighting makes reading a lot easier.

Second, you maybe don't need to define CohFrict-Burger interactions if burger inherits from CohFrict.
They will be calculated automatically as CohFrict-CohFrict. If it's not a good behavior then indeed you need this third mixed type, else you don't.

I paste it here. How can I attach files?

// 2014 © Behzad Majidi <email address hidden>
// 2014 © Ricardo Pieralisi <email address hidden>

#pragma once

#include<core/Material.hpp>
#include<pkg/dem/FrictPhys.hpp>
#include<pkg/dem/ViscoelasticPM.hpp>
#include<pkg/dem/CohesiveFrictionalContactLaw.hpp>
#include<pkg/common/Dispatching.hpp>
#include<pkg/dem/ScGeom.hpp>
#include<pkg/dem/DemXDofGeom.hpp>
#include<pkg/common/MatchMaker.hpp>

//********************** CohBurgersMat ****************************/

class CohBurgersMat : public CohFrictMat
{
 public :
  virtual ~CohBurgersMat () {};

/// Serialization
 YADE_CLASS_BASE_DOC_ATTRS_CTOR(CohBurgersMat,CohFrictMat,"",
  ((Real,kmn,NaN,,"Normal elastic stiffness for Maxwell."))
  ((Real,cmn,NaN,,"Normal viscous constant for Maxwell."))
  ((Real,kms,NaN,,"Shear elastic stiffness for Maxwell."))
  ((Real,cms,NaN,,"Shear viscous constant for Maxwell."))
  ((Real,kkn,NaN,,"Normal elastic stiffness for Kelvin."))
  ((Real,ckn,NaN,,"Normal viscous constant for Kelvin."))
  ((Real,kks,NaN,,"Shear elastic stiffness for Kelvin."))
  ((Real,cks,NaN,,"Shear viscous constant for Kelvin.")),
  createIndex();
  );
/// Indexable
 REGISTER_CLASS_INDEX(CohBurgersMat,CohFrictMat);
};

REGISTER_SERIALIZABLE(CohBurgersMat);

//********************** CohBurgersPhys ****************************/

class CohBurgersPhys : public CohFrictPhys
{
 public :
  virtual ~CohBurgersPhys() {};
  void SetBreakingState() {cohesionBroken = true; normalAdhesion = 0; shearAdhesion = 0;};

 YADE_CLASS_BASE_DOC_ATTRS_CTOR(CohBurgersPhys,CohFrictPhys,"",
  ((Real,kmn,NaN,,"Stiffness of Maxwell's spring(normal)"))
  ((Real,kkn,NaN,,"Stiffness of Kelvin's spring(normal)"))
  ((Real,cmn,NaN,,"Viscosity of Maxwell's dashpot(normal)"))
      ((Real,ckn,NaN,,"Viscosity of Kelvin's dashpot(normal)"))
  ((Real,kms,NaN,,"Stiffness of Maxwell's spring(shear)"))
  ((Real,kks,NaN,,"Stiffness of Kelvin's spring(shear)"))
  ((Real,cms,NaN,,"Viscosity of Maxwell's dashpot(shear)"))
      ((Real,cks,NaN,,"Viscosity of Kelvin's dashpot(shear)"))
  ((Vector3r,ukn,Vector3r(0,0,0),,"Normal displacement"))
  ((Vector3r,uks,Vector3r(0,0,0),,"Shear displacement")),
  createIndex();
 )
 REGISTER_CLASS_INDEX(CohBurgersPhys,CohFrictPhys);
};

REGISTER_SERIALIZABLE(CohBurgersPhys);

//********************** Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys ****************************/

class Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys: public IPhysFunctor
 {
 public :
  virtual void go( const shared_ptr<Material>& b1,
      const shared_ptr<Material>& b2,
      const shared_ptr<Interaction>& interaction);
  virtual void Calculate_CohBurgersMat_CohBurgersMat_CohBurgersPhys(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction, shared_ptr<CohBurgersPhys> phys);
  int cohesionDefinitionIteration;
 YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys,IPhysFunctor,
 "Convert 2 instances of :yref:`CohBurgersMat` to :yref:`CohBurgersPhys` using the rule of consecutive connection.",
 ((bool,setCohesionNow,false,,"If true, assign cohesion to all existing contacts"))
 ((bool,setCohesionOnNewContacts,false,,"If true, assign cohesion at all new contacts."))
 ,
 cohesionDefinitionIteration = -1;
 );
 FUNCTOR2D(CohBurgersMat,CohBurgersMat);
};
REGISTER_SERIALIZABLE(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys);

class Ip2_FrictMat_CohBurgersMat_CohBurgersPhys: public IPhysFunctor
{
 public :
  virtual void go( const shared_ptr<Material>& b1,
     const shared_ptr<Material>& b2,
     const shared_ptr<Interaction>& interaction);
  virtual void Calculate_FrictMat_CohBurgersMat_CohBurgersPhys(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction, shared_ptr<CohBurgersPhys> phys);
  int cohesionDefinitionIteration;
 YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_FrictMat_CohBurgersMat_CohBurgersPhys,IPhysFunctor,
 "converts an interaction of FrictMat and CohBurgersMat to CohBurgersPhys",
 ((bool,setCohesionNow,false,,"If true, assign cohesion to all existing contacts"))
 ((bool,setCohesionOnNewContacts,false,,"If true, assign cohesion at all new contacts."))
 ,
 cohesionDefinitionIteration = -1;
 );
 FUNCTOR2D(FrictMat,CohBurgersMat);
};

REGISTER_SERIALIZABLE(Ip2_FrictMat_CohBurgersMat_CohBurgersPhys);

//******************** Law2_ScGeom_CohBurgersPhys_CohesiveBurgers *************************/

class Law2_ScGeom_CohBurgersPhys_CohesiveBurgers: public LawFunctor {
 public:
  virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction*);

 FUNCTOR2D(ScGeom,CohBurgersPhys);

 YADE_CLASS_BASE_DOC(Law2_ScGeom_CohBurgersPhys_CohesiveBurgers, LawFunctor,"");

};
REGISTER_SERIALIZABLE(Law2_ScGeom_CohBurgersPhys_CohesiveBurgers);

bool computeForceTorqueCohBurgers(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I, Vector3r & force, Vector3r & torque1, Vector3r & torque2);

Could you also show us the Python script you use?
Thanks
Jan

2015-02-05 19:21 GMT+01:00 behzad <email address hidden>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> behzad posted a new comment:
> I paste it here. How can I attach files?
>
> // 2014 © Behzad Majidi <email address hidden>
> // 2014 © Ricardo Pieralisi <email address hidden>
>
> #pragma once
>
> #include<core/Material.hpp>
> #include<pkg/dem/FrictPhys.hpp>
> #include<pkg/dem/ViscoelasticPM.hpp>
> #include<pkg/dem/CohesiveFrictionalContactLaw.hpp>
> #include<pkg/common/Dispatching.hpp>
> #include<pkg/dem/ScGeom.hpp>
> #include<pkg/dem/DemXDofGeom.hpp>
> #include<pkg/common/MatchMaker.hpp>
>
>
> //********************** CohBurgersMat ****************************/
>
> class CohBurgersMat : public CohFrictMat
> {
> public :
> virtual ~CohBurgersMat () {};
>
> /// Serialization
> YADE_CLASS_BASE_DOC_ATTRS_CTOR(CohBurgersMat,CohFrictMat,"",
> ((Real,kmn,NaN,,"Normal elastic stiffness for Maxwell."))
> ((Real,cmn,NaN,,"Normal viscous constant for Maxwell."))
> ((Real,kms,NaN,,"Shear elastic stiffness for Maxwell."))
> ((Real,cms,NaN,,"Shear viscous constant for Maxwell."))
> ((Real,kkn,NaN,,"Normal elastic stiffness for Kelvin."))
> ((Real,ckn,NaN,,"Normal viscous constant for Kelvin."))
> ((Real,kks,NaN,,"Shear elastic stiffness for Kelvin."))
> ((Real,cks,NaN,,"Shear viscous constant for Kelvin.")),
> createIndex();
> );
> /// Indexable
> REGISTER_CLASS_INDEX(CohBurgersMat,CohFrictMat);
> };
>
> REGISTER_SERIALIZABLE(CohBurgersMat);
>
>
> //********************** CohBurgersPhys ****************************/
>
> class CohBurgersPhys : public CohFrictPhys
> {
> public :
> virtual ~CohBurgersPhys() {};
> void SetBreakingState() {cohesionBroken = true;
> normalAdhesion = 0; shearAdhesion = 0;};
>
> YADE_CLASS_BASE_DOC_ATTRS_CTOR(CohBurgersPhys,CohFrictPhys,"",
> ((Real,kmn,NaN,,"Stiffness of Maxwell's spring(normal)"))
> ((Real,kkn,NaN,,"Stiffness of Kelvin's spring(normal)"))
> ((Real,cmn,NaN,,"Viscosity of Maxwell's dashpot(normal)"))
> ((Real,ckn,NaN,,"Viscosity of Kelvin's dashpot(normal)"))
> ((Real,kms,NaN,,"Stiffness of Maxwell's spring(shear)"))
> ((Real,kks,NaN,,"Stiffness of Kelvin's spring(shear)"))
> ((Real,cms,NaN,,"Viscosity of Maxwell's dashpot(shear)"))
> ((Real,cks,NaN,,"Viscosity of Kelvin's dashpot(shear)"))
> ((Vector3r,ukn,Vector3r(0,0,0),,"Normal displacement"))
> ((Vector3r,uks,Vector3r(0,0,0),,"Shear displacement")),
> createIndex();
> )
> REGISTER_CLASS_INDEX(CohBurgersPhys,CohFrictPhys);
> };
>
> REGISTER_SERIALIZABLE(CohBurgersPhys);
>
>
> //********************** Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys
> ****************************/
>
> class Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys: public IPhysFunctor
> {
> public :
> virtual void go( const
> shared_ptr<Material>& b1,
> const
> shared_ptr<Material>& b2,
> const
> shared_ptr<Interaction>& interaction);
> virtual void
> Calculate_CohBurgersMat_CohBurgersMat_CohBurgersPhys(const
> shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const
> shared_ptr<Interaction>& interaction, shared_ptr<CohBurgersPhys> phys);
> int cohesionDefinitionIteration;
>
> YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys,IPhysFunctor,
> "Convert 2 instances of :yref:`CohBurgersMat` to
> :yref:`CohBurgersPhys` using the rule of consecutive connection.",
> ((bool,setCohesionNow,false,,"If true, assign cohesion to all
> existing contacts"))
> ((bool,setCohesionOnNewContacts,false,,"If true, assign cohesion
> at all new contacts."))
> ,
> cohesionDefinitionIteration = -1;
> );
> FUNCTOR2D(CohBurgersMat,CohBurgersMat);
> };
> REGISTER_SERIALIZABLE(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys);
>
>
>
> class Ip2_FrictMat_CohBurgersMat_CohBurgersPhys: public IPhysFunctor
> {
> public :
> virtual void go( const shared_ptr<Material>& b1,
> const shared_ptr<Material>& b2,
> const shared_ptr<Interaction>&
> interaction);
> virtual void
> Calculate_FrictMat_CohBurgersMat_CohBurgersPhys(const shared_ptr<Material>&
> b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>&
> interaction, shared_ptr<CohBurgersPhys> phys);
> int cohesionDefinitionIteration;
>
> YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_FrictMat_CohBurgersMat_CohBurgersPhys,IPhysFunctor,
> "converts an interaction of FrictMat and CohBurgersMat to
> CohBurgersPhys",
> ((bool,setCohesionNow,false,,"If true, assign cohesion to all
> existing contacts"))
> ((bool,setCohesionOnNewContacts,false,,"If true, assign cohesion
> at all new contacts."))
> ,
> cohesionDefinitionIteration = -1;
> );
> FUNCTOR2D(FrictMat,CohBurgersMat);
> };
>
> REGISTER_SERIALIZABLE(Ip2_FrictMat_CohBurgersMat_CohBurgersPhys);
>
>
> //******************** Law2_ScGeom_CohBurgersPhys_CohesiveBurgers
> *************************/
>
>
> class Law2_ScGeom_CohBurgersPhys_CohesiveBurgers: public LawFunctor {
> public:
> virtual bool go(shared_ptr<IGeom>& _geom,
> shared_ptr<IPhys>& _phys, Interaction*);
>
> FUNCTOR2D(ScGeom,CohBurgersPhys);
>
> YADE_CLASS_BASE_DOC(Law2_ScGeom_CohBurgersPhys_CohesiveBurgers,
> LawFunctor,"");
>
> };
> REGISTER_SERIALIZABLE(Law2_ScGeom_CohBurgersPhys_CohesiveBurgers);
>
> bool computeForceTorqueCohBurgers(shared_ptr<IGeom>& _geom,
> shared_ptr<IPhys>& _phys, Interaction* I, Vector3r & force, Vector3r &
> torque1, Vector3r & torque2);
>
> --
> 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
>

Sure!

Here it is:

The first one gives CohBurgers-CohBurgers interaction and the model prediction coincides the analytical response.
However, the problem is the FrictMat-CohBurgersMat interaction. (the second python script)

========CohBurgers-CohBurgers interaction===============================

O.reset()
from yade import utils, plot
from yade import pack, qt

id_Mat1=O.materials.append(CohBurgersMat(kmn=1.6e8,kkn=1.2e8,cmn=3e8,ckn=3e8,
kms=1.2e8,kks=1.2e8,cms=3e8,cks=3e8,normalCohesion= 1e5, shearCohesion= 1e3, isCohesive= True, young=1e2,
density=2600, poisson=0.3, frictionAngle= 0.4))
Mat=O.materials[id_Mat1]

s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat)

O.bodies.append(s1)
O.bodies.append(s2)

O.engines=[
ForceResetter(),
ForceEngine(force=(0,0,7.0e4),ids=[1],label='fEngine'),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
[Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys(setCohesionNow=True)],
[Law2_ScGeom_CohBurgersPhys_CohesiveBurgers(),Law2_ScGeom_FrictPhys_CundallStrack()]
),
NewtonIntegrator(damping=0.7,gravity=[0,0,0]),
PyRunner(command="fEngine.force=(0,0,0)", iterPeriod=10000000),
PyRunner(command='AutoData()',iterPeriod=100000)
]

def AutoData():
 f_app=7.0e4
 t_loading=10
 cmn=3e8
 kmn=1.6e8
 kkn=1.2e8
 ckn=3e8
 j1= 1/kmn
 j2= 1/kkn
 j3= kkn/ckn
 j4= 1/cmn
 a1= cmn/kmn+cmn*(1.0/kmn+1/kmn)
 a2=cmn*cmn/(kmn*kmn)
 b1= cmn
 b2=cmn*cmn/kmn
 z1=(-a1+sqrt(a1*a1-4.0*a2))/(2.0*a2)
 z2=(-a1-sqrt(a1*a1-4.0*a2))/(2.0*a2)
 aa1=(b2*z1+b1)/a2/(z1-z2)
 aa2=(b2*z2+b1)/a2/(z2-z1)
 if O.time < t_loading:
  displacement=f_app*(j1+j2*(1-exp(-(O.time)*j3))+(O.time)*j4)
 else:
  displacement=f_app*((t_loading/cmn)+j2*(exp(-(O.time)*j3))*(exp((t_loading)*j3)-1));

        plot.addData(t=O.time,d_numerical=0.02-(O.bodies[s2.id].state.pos[2]),d_analytical=-(displacement))

plot.plots={'t':(('d_numerical',('d_analytical','^:'),))}
plot.plot()

O.dt=1e-6

qt.View()

===================================================================

========FrictMat-CohBurgers Interaction================================

O.reset()
from yade import utils, plot
from yade import pack, qt

id_Mat1=O.materials.append(CohBurgersMat(kmn=1.6e8,kkn=1.2e8,cmn=3e8,ckn=3e8,
kms=1.2e8,kks=1.2e8,cms=3e8,cks=3e8,normalCohesion= 1e3, shearCohesion= 1e3, isCohesive= True, young=1e2,
density=2600, poisson=0.3, frictionAngle= 0.4))
Mat1=O.materials[id_Mat1]

id_Mat2=O.materials.append(FrictMat(young=1e8,poisson=0.3,density=2600,frictionAngle=1))
Mat2=O.materials[id_Mat2]

s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)

O.bodies.append(s1)
O.bodies.append(s2)

O.engines=[
ForceResetter(),
ForceEngine(force=(0,0,-2e-4),ids=[1],label='fEngine'),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
[Ip2_FrictMat_CohBurgersMat_CohBurgersPhys()],
[Law2_ScGeom_CohBurgersPhys_CohesiveBurgers(),Law2_ScGeom_FrictPhys_CundallStrack()]
),
NewtonIntegrator(damping=0.7,gravity=[0,0,0]),
PyRunner(command="fEngine.force=(0,0,0)", iterPeriod=5000000),
PyRunner(command='AutoData()',iterPeriod=150000)
]

def AutoData():
 f_app=2e-4
 t_loading=5
 cmn=3e8
 kmn=0.615384615e8
 kkn=1.2e8
 ckn=3e8
 j1= 1/kmn
 j2= 1/kkn
 j3= kkn/ckn
 j4= 1/cmn
 a1= cmn/kmn+cmn*(1.0/kmn+1/kmn)
 a2=cmn*cmn/(kmn*kmn)
 b1= cmn
 b2=cmn*cmn/kmn
 z1=(-a1+sqrt(a1*a1-4.0*a2))/(2.0*a2)
 z2=(-a1-sqrt(a1*a1-4.0*a2))/(2.0*a2)
 aa1=(b2*z1+b1)/a2/(z1-z2)
 aa2=(b2*z2+b1)/a2/(z2-z1)
 if O.time < t_loading:
  displacement=f_app*(j1+j2*(1-exp(-(O.time)*j3))+(O.time)*j4)
 else:
  displacement=f_app*((t_loading/cmn)+j2*(exp(-(O.time)*j3))*(exp((t_loading)*j3)-1));

        plot.addData(t=O.time,d_numerical=0.02-(O.bodies[s2.id].state.pos[2]),d_analytical=displacement)

plot.plots={'t':(('d_numerical',('d_analytical','^:'),))}
plot.plot()

O.dt=1e-6

qt.View()

=================================================================

If this "kmn" attribute is the expected value at one time, and if your interaction does not have any physics at the beginning (your AttributeError message), I would say that your Ip2 c++ code does - at least, partly - correctly its job, but that maybe there is problems in other parts of the InteractionLoop..
Your python script would be useful, too, surely..

As a very basical advice, you might try to ouptut messages from the c++ code to check when the physics of the interaction is created, and if this time corresponds to what you think.

Hi Behzad,

I would continue by repeating what is actually the problem :-)
- interaction physics creation - at the beginning of the simulation, the
spheres are just touching, it is the matter of ronging errors if
interaction (and the physics too) are created or not
- difference of numerical and analytical results - no idea, can you check
if all parameters of O.interactions[0,1].phys are correct in Python?

Jan

2015-02-05 19:56 GMT+01:00 Jérôme Duriez <
<email address hidden>>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> Jérôme Duriez proposed the following answer:
> If this "kmn" attribute is the expected value at one time, and if your
> interaction does not have any physics at the beginning (your AttributeError
> message), I would say that your Ip2 c++ code does - at least, partly -
> correctly its job, but that maybe there is problems in other parts of the
> InteractionLoop..
> Your python script would be useful, too, surely..
>
> As a very basical advice, you might try to ouptut messages from the c++
> code to check when the physics of the interaction is created, and if
> this time corresponds to what you think.
>
> --
> 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 see that "Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys" and "Ip2_FrictMat_CohBurgersMat_CohBurgersPhys" are almost identical. There is only a difference in how Kmn is calculated when FrictMat is in contact with a BurgersMat. All other calculations and statements are the same. So, I'm wondering where the problem is for FrictMat-CohBurgersMat contact.

Something which can be of problem is how the materials sequence can affect the interaction here and I don't know if I've written the code correctly. I'm talking about:

const shared_ptr<FrictMat>& mat1 = YADE_PTR_CAST<FrictMat>(b1);
const shared_ptr<CohBurgersMat>& mat2 = YADE_PTR_CAST<CohBurgersMat>(b2);

in void Ip2_FrictMat_CohBurgersMat_CohBurgersPhys::Calculate_FrictMat_CohBurgersMat_ of cpp file.

Is this a correct statement? I mean, if in my Python file we switch the materials of O.bodies[0] and O.bodies[1], does it affect the interaction? It seems it does! I did this and by running the code, O.bodies[1] disappears and O.bodies[1].state.pos gives:
 -> [4]: Vector3(nan,nan,nan)

Hi Behzad,
the non-spherical clumps should be possible to use in latest trunk version
(r3576). I added an example, in principle doing exactly what you want,
yade/examples/test/clump-facet.py
cheers
Jan

2015-02-05 20:11 GMT+01:00 Jan Stránský <
<email address hidden>>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> Jan Stránský proposed the following answer:
> Hi Behzad,
>
> I would continue by repeating what is actually the problem :-)
> - interaction physics creation - at the beginning of the simulation, the
> spheres are just touching, it is the matter of ronging errors if
> interaction (and the physics too) are created or not
> - difference of numerical and analytical results - no idea, can you check
> if all parameters of O.interactions[0,1].phys are correct in Python?
>
> Jan
>
>
>
> 2015-02-05 19:56 GMT+01:00 Jérôme Duriez <
> <email address hidden>>:
>
> > Question #261724 on Yade changed:
> > https://answers.launchpad.net/yade/+question/261724
> >
> > Jérôme Duriez proposed the following answer:
> > If this "kmn" attribute is the expected value at one time, and if your
> > interaction does not have any physics at the beginning (your
> AttributeError
> > message), I would say that your Ip2 c++ code does - at least, partly -
> > correctly its job, but that maybe there is problems in other parts of the
> > InteractionLoop..
> > Your python script would be useful, too, surely..
> >
> > As a very basical advice, you might try to ouptut messages from the c++
> > code to check when the physics of the interaction is created, and if
> > this time corresponds to what you think.
> >
> > --
> > 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 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, wrong thread.. please ignore my last message
Jan

2015-02-07 19:01 GMT+01:00 Jan Stránský <
<email address hidden>>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> Jan Stránský proposed the following answer:
> Hi Behzad,
> the non-spherical clumps should be possible to use in latest trunk version
> (r3576). I added an example, in principle doing exactly what you want,
> yade/examples/test/clump-facet.py
> cheers
> Jan
>
>
> 2015-02-05 20:11 GMT+01:00 Jan Stránský <
> <email address hidden>>:
>
> > Question #261724 on Yade changed:
> > https://answers.launchpad.net/yade/+question/261724
> >
> > Jan Stránský proposed the following answer:
> > Hi Behzad,
> >
> > I would continue by repeating what is actually the problem :-)
> > - interaction physics creation - at the beginning of the simulation, the
> > spheres are just touching, it is the matter of ronging errors if
> > interaction (and the physics too) are created or not
> > - difference of numerical and analytical results - no idea, can you check
> > if all parameters of O.interactions[0,1].phys are correct in Python?
> >
> > Jan
> >
> >
> >
> > 2015-02-05 19:56 GMT+01:00 Jérôme Duriez <
> > <email address hidden>>:
> >
> > > Question #261724 on Yade changed:
> > > https://answers.launchpad.net/yade/+question/261724
> > >
> > > Jérôme Duriez proposed the following answer:
> > > If this "kmn" attribute is the expected value at one time, and if your
> > > interaction does not have any physics at the beginning (your
> > AttributeError
> > > message), I would say that your Ip2 c++ code does - at least, partly -
> > > correctly its job, but that maybe there is problems in other parts of
> the
> > > InteractionLoop..
> > > Your python script would be useful, too, surely..
> > >
> > > As a very basical advice, you might try to ouptut messages from the c++
> > > code to check when the physics of the interaction is created, and if
> > > this time corresponds to what you think.
> > >
> > > --
> > > 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 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 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
>

Problem is in Ip2_FrictMat_CohBurgersMat_CohBurgersPhys.

The problem is in my model.py, when the 1st body is FrictMat and the 2nd body is CohBurgersMat everything is working well. But, if it's the reverse (1st body is CohBurgersMat and the 2nd body is FrictMat ) the interaction is Nan.

This means in Ip2 materials are not well defined. But, I don't know how can I define materials so the interaction is formed whenever a FrictMat is in contact with a CohBurgersMat whatever the sequence is.

Regards,
Behzad

I think you need to define Ip2::goReverse for dealing with the failing sequence.
Grep "goReverse" in the source code, you will find examples of that.
Bruno

Goo idea

However, I tried it the same way it appears in the examples, but it didn't work!

I have go function in cpp as:

void Ip2_FrictMat_CohBurgersMat_CohBurgersPhys::go(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction)
{

  if(interaction->phys) return;
  shared_ptr<CohBurgersPhys> phys (new CohBurgersPhys());

// const shared_ptr<FrictMat>& mat1 = YADE_PTR_CAST<FrictMat>(b1);
// const shared_ptr<CohBurgersMat>& mat2 = YADE_PTR_CAST<CohBurgersMat>(b2);

  FrictMat* mat1 = static_cast<FrictMat*>(b1.get());
  CohBurgersMat* mat2 = static_cast<CohBurgersMat*>(b2.get());

  ScGeom* geom= YADE_CAST<ScGeom*>(interaction->geom.get());

  Real mass1 = 1.0;
  Real mass2 = 1.0;
  Real Ra = geom->radius1;
  Real Rb = geom->radius2;

  mass1=Body::byId(interaction->getId1())->state->mass;
  mass2=Body::byId(interaction->getId2())->state->mass;

//Load the rheological parameters

  Real kn1 = (mat1->young); Real ks1 = (mat1->young)*(mat1->poisson);
  Real kmn2 = (mat2->kmn); Real kkn2 = (mat2->kkn);
  Real cmn2 = (mat2->cmn); Real ckn2 = (mat2->ckn);
  Real kms2 = (mat2->kms); Real kks2 = (mat2->kks);
  Real cms2 = (mat2->cms); Real cks2 = (mat2->cks);

  phys->kmn = contactParameterCalculation(kn1,kmn2);
  phys->kms = contactParameterCalculation(ks1,kms2);
  phys->cmn = cmn2;
  phys->cms = cms2;

  phys->kkn = kkn2;
  phys->kks = kks2;
  phys->ckn = ckn2;
  phys->cks = cks2;
  phys->MaxDef = 0.25*(Ra+Rb);
  phys->normalAdhesion= mat2->normalCohesion*0.5*pow(std::min(Ra,Rb),2);

  interaction->phys = shared_ptr<CohBurgersPhys>(phys);

  phys->cohesionBroken = false;
  phys->initCohesion = true;
  cohesionDefinitionIteration = -1;
  setCohesionNow = true;

}

and th goReverse as:

void Ip2_FrictMat_CohBurgersMat_CohBurgersPhys::goReverse(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction)
{
  interaction->swapOrder();
  go(b1,b2,interaction);

}

Hey guys,

This is a serious issue for me now. Any idea is highly appreciated.
I just cannot make this FrictMat_ConBurgersMat interaction. I mean it only works if id1 is FrictMat and id2 is ConBurgersMat and the reverse is making interaction nan.

Cheers

Probably a mistake in goReverse. I can't check details now sorry.
Best way to debug is to not use static_cast's, they hide the problems.
Replace them with YADE_CAST, and of course make a debug build.

Bruno

If this is a really really bad issue, maybe handle Frict_Burger as Frict_Frict in the short term. Or just don't mix Burger and frictional materials.
Isn't it possible to reproduce the Frict behaviour with the Burger classes?!
B

- Yes, I was able to handle it as Frict_Frict. But, it doesn't not really fit the physics of the contact in reality.

-replacing the static_cast with YADE_CAST didn't change anything.

- As for that last option, yeah, I will try to mimic FrictMat with Burgers class. But this goReverese thing seems to be an easy method and we have the examples in the source code also. I really don't know why it doesn't work here.

Thanks Bruno

>replacing the static_cast with YADE_CAST didn't change anything.

Just a good practice when trying to solve such problem. You are using a debug build, yes?

I tried:

cmake -DINSTALL_PREFIX=~/yade-trunk/ ~/yade-trunk/trunk/ -DCMAKE_BUILD_TYPE=DEBUG

As for myself, I did not know the "CMAKE_BUILD_TYPE" option... I am using the "DEBUG" one ("-DDEBUG=ON", then "make install"...)

Where did you find the option you're using ? Not on https://yade-dem.org/doc/installation.html#compilation, isn't it ?

No it wasn't on www.yade-dem.org.
I found it on a web page for make compilation.

Hi Behzad, could you place your complete .hpp and .cpp files together with
a python (preferably short) script reproducing the error somewhere on the
internet such that we can try to compile and test yade with your extension?
Thanks
Jan

2015-02-20 19:21 GMT+01:00 behzad <email address hidden>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> behzad posted a new comment:
> No it wasn't on www.yade-dem.org.
> I found it on a web page for make compilation.
>
> --
> 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
>

Hi Jan,

Here you go!

https://www.dropbox.com/s/xs5jr936d2lamus/CohBurgers.7z?dl=0

You've got the hpp and cpp files and a script file showing an example of interactions of CohBurgersMat and also a readme file which explains the material and the problem.

Cheers,
Behzad

Hi Behzad,

thanks for files. I compiled them and used your py file, with no problem at
all :-)

cheers
Jan

######################################################################
O.reset()
from yade import utils, plot
from yade import pack, qt

id_Mat1=O.materials.append(CohBurgersMat(kmn=2.0e8,kkn=1.0e8,cmn=3e8,ckn=3e8,
kms=2.0e8,kks=1.0e8,cms=3e8,cks=3e8,normalCohesion= 1e13, shearCohesion=
1e13, isCohesive= True, young=1e2,
density=2600, poisson=0.3, frictionAngle= 0.4))
Mat1=O.materials[id_Mat1]

id_Mat2=O.materials.append(FrictMat(young=1.0e8,poisson=0.3,density=1000,frictionAngle=1))
Mat2=O.materials[id_Mat2]

if 0: # works for both 0 and 1
s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)
else:
s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat1)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat2)

O.bodies.append(s1)
O.bodies.append(s2)

O.engines=[
ForceResetter(),
ForceEngine(force=(0,0,-2e4),ids=[1],label='fEngine'),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
[Ip2_FrictMat_CohBurgersMat_CohBurgersPhys(),Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys()],
[Law2_ScGeom_CohBurgersPhys_CohesiveBurgers()]
),
NewtonIntegrator(damping=0.7,gravity=[0,0,0]),
PyRunner(command="fEngine.force=(0,0,0)", iterPeriod=5000000),
PyRunner(command='AutoData()',iterPeriod=150000)
]

def AutoData():
f_app=2e4
t_loading=5
cmn=3e8
# kmn=2.0e8
kmn=0.66666667e8
kkn=1.0e8
ckn=3e8
j1= 1/kmn
j2= 1/kkn
j3= kkn/ckn
j4= 1/cmn
a1= cmn/kmn+cmn*(1.0/kmn+1/kmn)
a2=cmn*cmn/(kmn*kmn)
b1= cmn
b2=cmn*cmn/kmn
z1=(-a1+sqrt(a1*a1-4.0*a2))/(2.0*a2)
z2=(-a1-sqrt(a1*a1-4.0*a2))/(2.0*a2)
aa1=(b2*z1+b1)/a2/(z1-z2)
aa2=(b2*z2+b1)/a2/(z2-z1)
if O.time < t_loading:
displacement=f_app*(j1+j2*(1-exp(-(O.time)*j3))+(O.time)*j4)
else:
displacement=f_app*((t_loading/cmn)+j2*(exp(-(O.time)*j3))*(exp((t_loading)*j3)-1));
        plot.addData(t=O.time,d_numerical=0.02-(O.bodies[s2.id
].state.pos[2]),d_analytical=displacement)

plot.plots={'t':(('d_numerical',('d_analytical','^:'),))}
plot.plot()

O.dt=1e-6

qt.View()
######################################################################

2015-02-20 21:26 GMT+01:00 behzad <email address hidden>:

> Question #261724 on Yade changed:
> https://answers.launchpad.net/yade/+question/261724
>
> behzad posted a new comment:
> Hi Jan,
>
> Here you go!
>
> https://www.dropbox.com/s/xs5jr936d2lamus/CohBurgers.7z?dl=0
>
> You've got the hpp and cpp files and a script file showing an example of
> interactions of CohBurgersMat and also a readme file which explains the
> material and the problem.
>
> Cheers,
> Behzad
>
> --
> 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
>

Look, I don't get what you're doing by:

if 0: # works for both 0 and 1
s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)
else:
s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat1)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat2)

But if you run the file with

s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)

it works!
And if you change it to:

s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat1)
s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat2)

It does not!

>
>
> Look, I don't get what you're doing by:
>
>
> if 0: # works for both 0 and 1
> s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
> s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)
> else:
> s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat1)
> s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat2)
>
>
> But if you run the file with
>
> s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat2)
> s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat1)
>
>
> it works!
> And if you change it to:
>
> s1=utils.sphere([0.0,0.0,0.0],0.01,fixed=True,material=Mat1)
> s2=utils.sphere([0.0,0.0,2.0e-2],0.01,fixed=False,material=Mat2)
>
> It does not!
>

my if else is just putting your two version of s1= s2= to one file for easy
change of them

sorry, it does not work for your original files, I tried some modifications
and it was the effect of undo/redo that I had the impression.

I you modify your .cpp file, specifically

  FrictMat* mat1 = YADE_CAST<FrictMat*>(b1.get());
  CohBurgersMat* mat2 = YADE_CAST<CohBurgersMat*>(b2.get());

to

  int i1 = b1->getClassIndex(); // get actual material index of b1
  int i2 = b2->getClassIndex(); // get actual material index of b2
  int cbmi = CohBurgersMat::getClassIndexStatic(); // get index of
CohBurgersMat
  FrictMat* mat1;
  CohBurgersMat* mat2;
  if (i2 == cbmi) { // b2 is CohburgersMat and b1 is FrictMat
  mat1 = dynamic_cast<FrictMat*>(b1.get());
  mat2 = dynamic_cast<CohBurgersMat*>(b2.get());
  } else if (i1 == cbmi) { // b1 is CohburgersMat and b2 is FrictMat
  mat1 = dynamic_cast<FrictMat*>(b2.get());
  mat2 = dynamic_cast<CohBurgersMat*>(b1.get());
  } else { // should not happen, but to be sure..
  LOG_FATAL("TODO");
  }

it should work. If not, let me know.
Basically you check the types of material instances and switch materials if
needed
Cheers
Jan

Great!
yes, this is what I needed to have in the code! Thank you, Jan. I appreciate it!

cheers,
Behzad