Fixed step size and end-of-step event

[Simone Hendriks]

Hi,

I have two questions about the integration of the EOM:

a) The system under study implements a force F=-A*sign(v), where A is a constant and v is the speed of the point where the force acts. The actual implementation smoothes the sign() function to avoid integration issues. The problem of the latter is that it introduces an additional parameter (a "smoothing stiffness") which is not always easy to set, and it leads anyway to small step sizes when the stiffness parameter is large. CPODES already performs better for this particular system, but my question is, would it make sense to use the unsmoothed sign() function and fix the integration step size? I am not particularly concerned about the sign() transition, but I need to avoid the integrator to reduce the step size. For the most part, this system has v around 0, which leads to too small step sizes for most of the integration. How would you approach this?

b) Unrelated to the previous question, is there a way to call a function after an integration step has completed? I am aware of the auto-update state variables, but for a particular application it would be easier to do it so.

Thank you in advance!

Best regards,
Simone

[Michael Sherman]

Hi, Simone.

Apologies for the delayed reply.

(b) There is currently no provision in Simbody to call a function at step completion. This would be a nice addition to the TimeStepper class in case anyone would like to submit a Simbody PR!

(a) This is a great question. It's not clear to me whether a fixed step would be stable -- worth a try. You can of course improve the integration speed by relaxing the drift parameter, at the expense of stiction that doesn't completely stick. That's the most common fix but might not work for you. One possibility (depending on the specifics of your application) would be to try Clay Anderson's recently added "ExponentialSpring" contact model. That uses extra state variables to deal with stiction in a non-stiff way that integrators like. It is limited in its applicability, though -- see the documentation in the header here. (You would need to build Simbody from source to get that.) I don't have any experience with that myself so don't know how effective that might be. Would be good to hear if you give it a try, though, and I'm sure Clay would be thrilled to hear from you.

Regards,
Sherm

[Simone Hendriks]

Hi Sherm,

Thank you very much for the reply, and apologies for my late reply as well.

(b) I'm not yet confident about extending Simbody, but please if anyone reading interested in this feature wants to submit a PR, it would certainly come in handy.

(a) Clay's approach is indeed very interesting! While the class is not directly applicable in my case (F=-A*sign(v) is actually a torque between two bodies, and not related to any contact station), the theory behind could apply. If I understood correctly, the frictional term consists of two parts: a linear spring-damper (the term that actually avoids drift) and a pure damper, both saturated to mu*N. The rest position of the spring is updated when full sliding takes place, which depends on a specified transition velocity. Aside from the fact that this formulation can avoid drift, I am still unsure about the efficiency improvement here, because there is still a transition velocity that for some applications needs to be very small, but I guess the addition of a fictive spring allows to consider greater transition speeds, thus making the system less stiff? I'd like to give it a try as time allows and post the results, it will also be interesting to compare it against a smoothed approach.

Thank you again for the valuable insights!

Best regards,
Simone