contact stiffness change y time

Asked by behzad

Hi,

I work on FrictPhys to make it behave like a viscoelastic contact. So, the idea is allow the kn and ks to change by the time of contact.

How do we include contact time in the contact physics?

let's say as an example we want

kn1=mat1->young*Da

phys->kn= kn1*contactTime

Thanks for the hints.

Cheers
Behzad

Question information

Language:
English Edit question
Status:
Answered
For:
Yade Edit question
Assignee:
No assignee Edit question
Last query:
Last reply:
Revision history for this message
Jan Stránský (honzik) said :
#1

Hi Behzad,

please be more specific on "time of contact" and "Da" and "contactTime"
variables in your example. Are they related to the O.iter of contact
creation, the time from currecnt iteration to iteration the contact was
created or something else?

Both cases I mentioned are easy to implement.

cheers
Jan

2016-02-12 23:22 GMT+01:00 behzad <email address hidden>:

> New question #285260 on Yade:
> https://answers.launchpad.net/yade/+question/285260
>
>
>
> Hi,
>
> I work on FrictPhys to make it behave like a viscoelastic contact. So, the
> idea is allow the kn and ks to change by the time of contact.
>
> How do we include contact time in the contact physics?
>
> let's say as an example we want
>
> kn1=mat1->young*Da
>
> phys->kn= kn1*contactTime
>
> Thanks for the hints.
>
> Cheers
> Behzad
>
> --
> You received this question notification because your team yade-users 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
>

Revision history for this message
behzad (behzad-majidi) said :
#2

Jan,

time of the contact = contactTime = duration of contact.

let's say O.time= t1 when the contact is created. Then what I mean is to implement the following equation :

phys->kn= kn1*deltaT

There is condition, if the bodies are approaching each other:

deltaT=t1-O.time

I know deltaT becomes negative, the main equation will be more complex, this is only an example.

and if the bodies are getting away each other:

deltaT=O.time-t1

Regards,
Behzad

Revision history for this message
Jan Stránský (honzik) said :
#3

Hi Behzad,

in Ip2 functor, in function Ip2::go, you alwas have something like

interaction->phys = shared_ptr<SomePhys>(new SomePhys());

At this point, you have access to O.time, so you can easily set (SomePhys
class has to have timeOfCreation attribute)

interaction->phys->timeOfCreation = O.time; // O.time will not work, but
some c++ equivalent

In Law2, you can again access O.time and phys.timeOfCreation, so you can do
whatever with these two values..

Is this what you mean (from your description I am still not sure :-)?

cheers
Jan

2016-02-16 16:42 GMT+01:00 behzad <email address hidden>:

> Question #285260 on Yade changed:
> https://answers.launchpad.net/yade/+question/285260
>
> behzad posted a new comment:
> Jan,
>
>
> time of the contact = contactTime = duration of contact.
>
> let's say O.time= t1 when the contact is created. Then what I mean is to
> implement the following equation :
>
> phys->kn= kn1*deltaT
>
> There is condition, if the bodies are approaching each other:
>
> deltaT=t1-O.time
>
> I know deltaT becomes negative, the main equation will be more complex,
> this is only an example.
>
>
> and if the bodies are getting away each other:
>
> deltaT=O.time-t1
>
>
> Regards,
> Behzad
>
> --
> You received this question notification because your team yade-users 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
>

Revision history for this message
behzad (behzad-majidi) said :
#4

Yes, that's what I meant.

1- O.time will not work, but some c++ equivalent

  That's the problem. I need to know what's that equivalent.

2- What about the relative velocity? I need to include that in the equations as well.
The stiffness of the contact must increase as the penetrationDepth is increasing and the vice versa.

Revision history for this message
behzad (behzad-majidi) said :
#5

and we must have a dynamic process here.
AFAIK, when we have interaction.geom then using the Ip2, we calculate and assign the contact parameters, like Kn, and this is constant for that interaction.
But, here, even if there is interaction.geom and interaction.phys we again calculate and update the phys.kn

so the contact.phys is getting updated in each iteraction.

Regards,

Revision history for this message
Jérôme Duriez (jduriez) said :
#6

Regarding 1- of #4, "scene->time" should work
The preliminary definition of "scene" as "const shared_ptr<Scene>& scene=Omega::instance().getScene();" might be necessary depending on which class you're into.

Regarding #5, updating continuously phys.kn is possible and meaningful from the Law2..::go()

Revision history for this message
Jan Stránský (honzik) said :
#7

And concerning relative velocity (and slimilar), have a look at ScGeom
class [1]
cheers
Jan

[1] https://github.com/yade/trunk/blob/master/pkg/dem/ScGeom.hpp#L42

2016-02-16 22:47 GMT+01:00 Jérôme Duriez <
<email address hidden>>:

> Question #285260 on Yade changed:
> https://answers.launchpad.net/yade/+question/285260
>
> Jérôme Duriez proposed the following answer:
> Regarding 1- of #4, "scene->time" should work
> The preliminary definition of "scene" as "const shared_ptr<Scene>&
> scene=Omega::instance().getScene();" might be necessary depending on which
> class you're into.
>
> Regarding #5, updating continuously phys.kn is possible and meaningful
> from the Law2..::go()
>
> --
> You received this question notification because your team yade-users 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
>

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#8

Out of curiosity: if the force-displacement relation is non linear, as you seem to suggest by mentionning changes of kn, then why do you need kn at all?
There is nothing wrong in writing fn=a*un^2, for instance. And it does not require to update "a" for every calculation of fn.

Besides, this statement is not correct: "using the Ip2, we calculate and assign the contact parameters, like Kn, and this is constant for that interaction". The Ip2 functors are called at every iteration and for every interaction, hence it is not against any rule to compute a different kn all the time (though pointless IMO). What is true is that _some_ Ip2 functors return immediatly when the ip is already defined, but it is not mandatory.

Revision history for this message
behzad (behzad-majidi) said :
#9

Bruno:

The idea is to get a viscoelastic model embedded in an elastic model. In this case, I'd like to have Burger's model embedded in CohFrictPhys.

Similar works have been already published using other DEM platforms:

[1] A.C. Collop et. al. "On the use of discrete element modelling to simulate the viscoelastic deformation behaviour of an idealized asphalt mixture"

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#10

>The idea is to get a viscoelastic model embedded in an elastic model.

What is the advantage of not writing simply a viscoelastic equation? There are examples of viscoelasticity in the existing laws, they don't trick any stiffness.

Revision history for this message
behzad (behzad-majidi) said :
#11

The intended model is Burger's model. I have already derived CohBurgersMat from CohFrictMat and CohBurgersPhys from CohFrictPhys. This model seems to be working very well.

However, I just noticed that there's a problem in shear/rotational forces. This makes the model unstable. I worked on it but I couldn't find the problem. Then, I decided to work on CohFrictMat to get the Burger's embedded on it. So, that's the reason!

I can share my model files with you. Actually I already did. I wanted to make it public, so everybody can see and modify it. I wasn't able to commit my changes with my account.

Revision history for this message
Anton Gladky (gladky-anton) said :
#12

2016-03-02 21:47 GMT+01:00 behzad <email address hidden>:
> I can share my model files with you. Actually I already did. I wanted to
> make it public, so everybody can see and modify it. I wasn't able to
> commit my changes with my account.

Is it a technical problem? Can we help with it?

Anton

Revision history for this message
behzad (behzad-majidi) said :
#13

Hi Anton,

There're the model files:

CohBurgers.cpp and CohBurgers.hpp and a test file Burgers_test.py:

https://www.dropbox.com/sh/vkqpimkrjr56yic/AADDihKwYjo7MQ0MFAZzZac5a?dl=0

Would you please commit it?

Revision history for this message
Jan Stránský (honzik) said :
#14

Hello Behzad,

if the model is unstable, it is usually a problem of time step. What
happens if you decrease the time step?

With the "stiffness" approach, would you get different forces / moments
then in "viscoelastic" approach? If so, why not to change the computation
directly in the "viscoelastic" approach?

cheers
Jan

2016-03-02 21:47 GMT+01:00 behzad <email address hidden>:

> Question #285260 on Yade changed:
> https://answers.launchpad.net/yade/+question/285260
>
> behzad posted a new comment:
>
> The intended model is Burger's model. I have already derived CohBurgersMat
> from CohFrictMat and CohBurgersPhys from CohFrictPhys. This model seems to
> be working very well.
>
> However, I just noticed that there's a problem in shear/rotational
> forces. This makes the model unstable. I worked on it but I couldn't
> find the problem. Then, I decided to work on CohFrictMat to get the
> Burger's embedded on it. So, that's the reason!
>
> I can share my model files with you. Actually I already did. I wanted to
> make it public, so everybody can see and modify it. I wasn't able to
> commit my changes with my account.
>
> --
> You received this question notification because your team yade-users 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
>

Revision history for this message
behzad (behzad-majidi) said :
#15

No, AFAIK, it's not the time-step issue. Because, the rotational velocities go crazy even with very small time-steps.

I have not defined rotational stiffnesses in the viscoelastic model. and I really couldn't figure out what's wrong. Thus, I thought, making CohFrictMat to behave like a viscoelastic material is an option. This method has been a well-known approach in the literature. I mentioned a reference in my previous comment (#9)

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#16

I just checked the reference in #9. IMHO the approach is flawn and its implementation is a nightmare (section 5 as whole is weird on other points to).
Eq. 1 is a pseudo-stiffness derived in the context of monotonic loading (the force is imposed and kept constant for a duration t). As soon as the loading changes in time it is not possible to define "t" properly, actually it would need a convolution of eq. 1 with the loading history... not something you want to do.

If a viscoelastic relation is implemented based on fig.2 then it will give eq.1 automaticaly, as a particular case.
The CohFrict model is actually fairly close to that, it's just the "C_infinity" dashpot which is not there (it is there in Viscoelastic contact model, OTOH, but then the other dashpot is missing).

Bruno

Revision history for this message
behzad (behzad-majidi) said :
#17

- I know, it's not the best way to deal the problem!
- Perhaps, troubleshooting of my current implementation of CohBurgersMat is the best way. However, I was not able to figure out what's making shear/rotational velocities go crazy.

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#18

I understand your view. On my side I'm reluctant to spend time on fixing something when I think it cannot work. :)
Understanding the cause of your problem needs a clear vision of the governing equations, yet the only thing I can tell about your equations is that they do NOT describe a Burger model.

Revision history for this message
behzad (behzad-majidi) said :
#19

You're talking about equations in the reference article in #9. Right?

So, let's make my classic Burger's model implementation, in the shared files in #13, public. So, you and everybody can see, try and modify it.

cheers!

Can you help with this problem?

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

To post a message you must log in.