New 6DOF Spring2 constraint spring limiting issue. Bug?

[LEgregius]

Around line 788 of btGeneric6DofSpring2Constraint is the following code.

Code: Select all

                //limit stiffness (the spring should not be sampled faster that the quarter of its angular frequency)
                if( 0.25 < angularfreq * dt)
                {
                        ks = BT_ONE / dt / dt / btScalar(16.0) / m;
                }

I believe this is incorrect, but I still have problems even if it is corrected. My understanding of this line is that if the angular frequency times the dt (the full step time, not the substep ) is greater than 0.25, then it limits ks to be a the spring frequency that will be 0.25.

Angular freq is defined as:

Code: Select all

btScalar m = mA > mB ? mB : mA;
...
btScalar angularfreq = sqrt(ks / m);

The problem is solving for ks in the previous equations yields

Code: Select all

ks = (0.25 * 0.25) / dt / dt * m; 

0.25 squared is the same as 1/16, so

Code: Select all

ks = BT_ONE / dt / dt / btScalar(16.0) * m;

Note the *m not /m. This seems to yield the intended result, however, I would also argue that the value should be based on the substep, not the full step, Furthermore, since it caps the spring all the time, it will yield incorrect results, which is what was happening to me, since I have sprung mass vehicle. I had to increase the step time just to make the vehicle not collapse on the bump stops.
I understand that means the vehicle may not have been stable, but since it uses the full step time and not the substep, I think this is probably incorrect as it does run stably if I disable this code.

[Flix]

I know nothing about the internals of this new constraint (and of contraints in general), but, based on the equations you wrote, the correction:

Code: Select all

ks = BT_ONE / dt / dt / btScalar(16.0) * m;
// Or: ks = BT_ONE*m/(dt*dt*btScalar(16.0));

makes sense.

BTW: Incidentally, in this bullet3 issue, https://github.com/bulletphysics/bullet3/issues/345 the GaborPuhr made a "division by mass" instead of a "mass multiplication" when calculating "c". I wonder if there's some relations...

As far as the "substep" problem I can't say anything... Are you referring to the "full" Bullet timestep, or something related to constraints ? Are the "full" substeps intended to be calculated at a different timestep than the main Bullet timestep ? Is this related to using variable timestep only ? [However, as I've told you, I know nothing about how constraints are implemented].

P.S. As a side note: I've seen in btGeneric6DofSpring2Constraint.h/cpp all these methods:

Code: Select all

// MatrixToEulerXYZ from http://www.geometrictools.com/LibFoundation/Mathematics/Wm4Matrix3.inl.html
[...]

I've been using them (and they're "MatrixFromEuler" counterparts) since I started using Bullet.
Maybe they could be placed somewhere inside the Bullet LinearMath library to make them more accessible.

[IresomPCH]

I have some questions about the Generic6DofSpring2Constraint too.
The used bullet version is 2.83.4, and I'm testing a mass-spring-mass system to find the parameter unit.


The current parameters are:
MLCP Dantzig solver @ 1000Hz stepTime(interpole) + @60Hz stepTime(internal default)
RigidBody 1 (top) mass 1kg, fixed position.
RigidBody 2 (bottom) mass 1kg, only connected to the spring.
Spring, setup using Generic6DofSpring2Constraint, connect body1 and body2, allow to move only in the vertical direction.
movable range = 0.3;
stiffness = 500.0;// variable
damping = 0.1;
erp, cfm = 0.1;

However, no matter what the stiffness is, the vibrate frequency doesn't change much.
I've tested stiffness=500, 5000, 50000, 500000...50000000, and the vibrate frequency are about 2.77Hz ~ 3Hz accordingly.
What's the stiffness unit of the constraint?

Another phenomenon happens when the mass ratio between body1 and body2 is large.

ex: body1=50kg, body2=0.01kg, stiffness = 50000.0, damping = 0.1;
The spring seems to accumulate energy as the mass-spring-mass system is dropped to the ground.
And the system jump higher and higher...

that's doesn't make sense. It should at least maintain in the same level.
The system will be stable when the mass ratio is small enough.

ex:
body1=50kg and body2=0.05kg, the system will keep jump to the same height.
body1=50kg and body2=0.5kg, the jump height is decreasing and go to rest.

And ideas?


----------------------------------------------------------
Update:

Previously I was using stepSimulation(stepTime); // stepTime = 0.001s, internalTimeStep = 1/60s
I've tested using stepSimulation(stepTime, 0); // stepTime = 0.001s, no internalTimeStep
And the Generic6DofSpring2Constraint behavior is close to a real spring, but somehow unstable. Higher stiffness makes a stiffer spring.
l'm still wondering the metric unit of stiffness...
It looks like a step frequency(no internalTimeStep) of 200Hz or even higher is necessary to simulate a good spring characteristic...

[LEgregius]

IresomPCH, are you saying before or after the change I made? If you step into the code, you can see where it will limit the spring value. If that if triggers, the value of the spring drops way down because of the math error, but it would drop it to a consistent value regardless of your spring value assuming the spring value is greater than 0.25 of the angular frequency in a tick.

[IresomPCH]

@LEgregius

I'm using the latest 2.83.4, according to the release note:
"New btGeneric6DofSpring2Constraint, replacing the old one. Thanks to Gabor PUHR and Tamas Umenhoffer. See the example browser under API or examples/Dof6Spring2Setup"


If the "spring value limit" mechanism exists, is there a better way to setup a stiffer spring?
Or at least act as the parameters changes according to a real metric unit?

btGeneric6DofSpring2Constraint is the one acting like the real spring the most...
I've tried btGeneric6DofSpringConstraint, Slider with motor, SliderSpring...none of the results are good enough

[Erwin Coumans]

I re-opened the issue related to btGeneric6DofSpring2Constraint
https://github.com/bulletphysics/bullet3/issues/345
not sure when someone has time to look into it. Please provide reproduction cases using modified Bullet examples, such as Bullet/examples/Constraints/Dof6Spring2Setup.cpp and attach them there.

Best to use the default btSequentialImpulseConstraintSolver or btNNCGConstraintSolver and a small time step without sub steps, and possibly double precision build.

Code: Select all

double deltaTime = 0.001;
world->stepSimulation(deltaTime,0);

Note that using large mass ratios and collisions and MLCP Dantzig solver is likely a different issue, and better not to mix problems.

[Gabor Puhr]

Hi,

A lot of things have been mentioned.

First:
It is true that at the stiffness clamping the value should be multiplied with the mass. Sorry for the inconvenience. A patch will be created shortly.

The solver (and becasue of that the constraints) don't know anything about the "full step time". So it always uses the substeps. (but it is called timeStep inside (or fps sometimes))

LEgregius is right this check is for stability, but it has to be reconsidered if this hardcoding is really necessary or let the user tune it to the end (where it will blow up eventually)

For IresomPCH:
The ERP and CFM values alter the behavior of the Spring heavly. I recomend for you to use ERP>=0.9 (but less than 1).
I usually use 0.999, but in this case the simulation became very jumpy so other parameters have to be set correctly.
Also CFM should be low. 0 is usually not a problem.
Please try commenting out the cited part and please try it with other solvers also.

For the big mass ratio: it is a well know problem for physics simulations that big mass ratio is not handled well.
The usual solution is that you decrease the ratio.
Also you can try the NNCG solver that is design to handle big mass ratios. It is also faster than the Dantzig usually, but generally not as accurate. Still in this case it is probably better. (although 1:1000 is too much for any kind of solver)

Please note that a physics simulator can only aproximate the ideal solution. It is easy to see that if you increase the magnitued of the forces in the system you will need smaller timestep to have the same precision. Stiffer spring means bigger forces. (Also stricter constraints causes bigger "forces".)

If the physical accuracy is necessary you have to sacrifice a little computing capacity. 1000Hz is not that rare. (double precision should be considered in this case probably though)

For seting up a spring: What is usually used is the angular frequency of a spring. (omega0)
stiffnes = m * sqr(omega0)

( I hope I used the m well here )

Best Regards
Gabor

[IresomPCH]

Hello,
is this issue solved in the 2.83.5 version?

[Erwin Coumans]

Hello,
is this issue solved in the 2.83.5 version?

Yes, it is fixed in latest version (and 2.83.5), thanks to Gabor Puhr.

[TheIgors]

Helllo

I use 2.83.6.

//limit stiffness (the spring should not be sampled faster that the quarter of its angular frequency)
if(limot->m_springStiffnessLimited && 0.25 < angularfreq * dt)
{
ks = BT_ONE / dt / dt / btScalar(16.0) * m;
}
//avoid damping that would blow up the spring
if(limot->m_springDampingLimited && kd * dt > m)
{
kd = m / dt;
}

The prob is the corrected values are independent from current "error" (deviation from equilibrium). As a result if a force (or gravity) is big enough then the spring just stops to react on bigger stiffness/damping values. Thus I have to turn off the correction.

The setup values is also a pain. As I read above and here they should be squared and then multiplied by mass. But it's inconvenient to control and makes user unhappy. Example: a force = 100 and gravity -100, applied to a cube with the constraint at bottom. What stiffness/damping values should he enter to get a rational result? (like 45 degrees rotated cube). It's absolute not intuitive, requires a series of experiments and dependent on cube's height. How can I simplify the setup for user?

Thx