Adding a method in the yade class files associated to a new contact law

Asked by Rioual on 2020-05-21


This question follows the previous answers and advices of Jan (
I want to add a method in my yade classes of a contact law for different tests.

- If I refer to the syntax for previous classes ( Cpmphys in ConcretePM.hpp):
They use .def() with all the parameters then the declaration of staticmethod
This doesn't exactly look like what Jan proposes in #5 of ...

- Concerning the implementation of the method in the .cpp file, here is how it looks for me copying the syntax used for ConcretePM.cpp:


bool Law2_ScGeom_VsmPhys_Vsm::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
 int id1 = contact->getId1(), id2 = contact->getId2();

 ScGeom* geom = static_cast<ScGeom*>(ig.get());
 VsmPhys* phys = static_cast<VsmPhys*>(ip.get());
 Omega* om = static_cast<Omega*>(ip.get());

 if (geom->penetrationDepth < 0) {
  if (neverErase) {
   phys->shearForce = Vector3r::Zero();
   phys->normalForce = Vector3r::Zero();
  } else
   return false;

 Real& Rav = 2 * (geom->refR1 * geom->refR2)/(geom->refR1 + geom->refR2)
 Real& mun = 3 * phys->viscos *(3.14159*Rav**2)*(phys->alfa)**(1/2)
 Real& mut = mun/3

 Real& Epsilonpn = geom->incidentVel * geom->normal * 1/(geom->refR1 + geom->refR2 - geom->penetrationDepth)*geom->normal

 Real& Fwakai = phys->Gamma * Rav * (-11.53 * (alfa) + 13.58 * (alfa)**(1/2) + 0.40)

 phys->normalForce = (mun * Epsilonpn - Fwakai)* geom->normal;

 phys->shearForce = mut*geom->shearIncrement()/om->dt

 phys->alfa = alfa + om->dt *( (1-alfa)/(alfa)**(1/2)*(0.38/phys->viscos)*(mun * Epsilonpn - Fwakai)/(3.14159*Rav**2) + (3/4*(3.14159)**(1/2))**(1/3)*(phys->surftens/phys->viscos)*((3.14159*Rav**2)**(-3/2) )

so how should i modify this eventually to implement Jan's option of definition of the method in .cpp ?? :
Real MyPhys :: computeNewAlpha(...) {
   // here you put the computation, wither directly or possibly using other functions, depends..
   return newAlpha;

I hope this is clear enough....

Thanks for your input,
All the best


Question information

English Edit question
Yade Edit question
No assignee Edit question
Last query:
Last reply:
Jan Stránský (honzik) said : #1

> They use .def() with all the parameters then ... staticmethod

all the parameters is not important, I think it is there for python docstrings and default parameters, you can use just
.def("pythonMethodName", &VsmPhys::someMethod)

> They use .def() ... then the declaration of staticmethod

it is not declaration, just conversion to Python. The &VsmPhys::someMethod has to be declared in the class before, like "static Real funcG(...)" in CpmPhys

> so how should i modify this ...
> I hope this is clear enough

sorry, I did not get the point :-)

specifically for CpmPhys::funcG, have a look how it is:
- declared in .hpp file [1]
- exposed to python [2,3] (as mentioned before, py::arg stuff is not needed)
- defined in .cpp file [4]
- used in Law2::go [5]

is that what you were asking?



Rioual (francois-rioual-v) said : #2

Yes, jan, thanks, this helps me a lot.
But concerning the class computeNewAlfa, i have a further question:

How do i precise in the formulation (hpp and cpp files) that OldAlfa is equal to NewAlfa for the previous time step corresponding to a specific interaction (how is the value of NewAlfa stored between two time steps... ?)..

I hope i am clear,
Thank you for your précisions,



Jan Stránský (honzik) said : #3

> concerning the class computeNewAlfa

it is a method/function, not class :-)

> How do i precise ...

however you want / your model needs..

> phys->alfa = alfa + om->dt *( ... )
> How do i precise ... that OldAlfa is equal to NewAlfa for the previous time step

depnding on the formulation, e.g. you can do it like in the original post, just the "phys->alfa = alfa + om->dt *( ... )" is defined in computeNewAlfa function (which can be composed of several other functions, like computeTerm1(...) etc.)

bool Law2_ScGeom_VsmPhys_Vsm::go(...)
   VsmPhys::computeAnything(phys->alfa); // phys->alfa is "oldAlfa"
   phys->alfa = computeNewAlfa(phys->alfa,...); // phys->alfa argument is "oldAlfa", phys->alfa result is "newAlfa"

or if you want the code being more "self-explanatory", the same can be rewritten as:

bool Law2_ScGeom_VsmPhys_Vsm::go(...)
   Real oldAlfa = phys->alfa;
   Real newAlfa = VsmPhys::computeNewAlfa(oldAlfa,...);
   // now you can do computations with both oldAlfa and newAlfa if needed
   phys->alfa = newAlfa;

> corresponding to a specific interaction

being member of IPhys = corresponding to a specific interaction

> how is the value of NewAlfa stored between two time steps... ?)

as phys->alfa (??)
if you need values in two time steps (current and previous), you can store two values in IPhys


Can you help with this problem?

Provide an answer of your own, or ask Rioual for more information if necessary.

To post a message you must log in.