Reference files for chaotic models should be short enough

[casella]

The MSL contains some chaotic models, which are theoretically guaranteed to provide diverging trajectories over time, no matter how accurate the tool and the integration algorithm is. Of course there is nothing wrong in simulating chaotic behaviour, but it doesn't make any sense to provide it as a reference, because it is theoretically impossible for any other tool than the one which produced the reference trajectory to reproduce it with a given degree of accuracy, if the simulation interval is long enough.

This issue has already been discussed, and there are proposals to add a special annotation to use a shorter StopTime for cross-check testing.

Until this is standardized, I would suggest to implement this manually, by re-running those cases with a shorter StopTime than the one found in the experiment annotation, which is meant to demonstrate the chaotic behaviour.

So far, I identified three such models:

• Modelica.Electrical.Analog.Examples.ChuaCircuit: classical example of chaotic electronic oscillator system
• Modelica.Mechanics.MultiBody.Examples.Elementary.PointGravityWithPointMasses2: N > 2 body systems can have chaotic motion
• Modelica.Mechanics.MultiBody.Examples.Rotational3DEffects.BevelGear1D: I don't know if this has ever analyzed from this point of view, but the motion is a complex nonlinear superposition of rotations, as in Furuta's pendulum, so I guess this may also be chaotic.

Based on the results obtained with OpenModelica, I would suggest to set StopTime to the following values when generating the reference trajectories of these three cases:

• BevelGear: 0.4 s
• ChuaCircuit: 4000 s
• PointGravityWithPointMasses2: 2s

This ensures sufficient data for regression testing and for cross-validation of tools, without incurring in the issue caused by chaotic motion.

If there is consensus on this, i can do it myself, or possibly @beutlich could do it.

[christiankral]

I am in favor of applying the proposed change to the ChuaCircuit reference trajectories

[MartinOtter]

Makes all sense, please apply the changes

[casella]

@beutlich, if you have the whole infrastructure ready and you can do that I would appreciate it, otherwise please let me know and I can take care of this myself. Thanks!

[HansOlsson]

To me the criteria for which models to select isn't clear. Is it just models where OpenModelica's result doesn't match the reference or was there some actual analysis of chaoticness?

I don't see that Modelica.Mechanics.MultiBody.Examples.Elementary.PointGravityWithPointMasses2 is a chaotic N-body system as the gravity field is fixed; and it just simulates one complicated rigid body in two different ways.

It still somewhat chaotic and becomes unpredictable after about 6s, but the default simulation interval is 3s. The much simpler Modelica.Mechanics.MultiBody.Examples.Elementary.PointGravity becomes unpredictable after 20-30s.

However, any rotating model with non-fixed frequency has a tendency for the phase to change in simulations causing references to be incorrect (and if it there are more degrees of freedom everything breaks then). That happens in BevelGear, but also in e.g., Modelica.Mechanics.MultiBody.Examples.Loops.EngineV6

[casella]

To me the criteria for which models to select isn't clear. Is it just models where OpenModelica's result doesn't match the reference or was there some actual analysis of chaoticness?

@HansOlsson, of course there was some actual analysis. Here it is in full detail for your reference.

Chua's circuit is a classical example of chaotic motion, see, e.g., here.

I don't see that Modelica.Mechanics.MultiBody.Examples.Elementary.PointGravityWithPointMasses2 is a chaotic N-body system as the gravity field is fixed; and it just simulates one complicated rigid body in two different ways.

As Newton demonstrated quite some time ago, a point mass body in a point gravity field moves along a periodic elliptical trajectory, according to Kepler's laws of motion. This is of course not chaotic. It is fairly straightforward to prove that a sphere of uniform density is gravitationally equivalent to a point mass located at its center, using Gauss's theorem and some symmetry arguments, so, no chaos in that case as well.

However, as soon as you have a less trivial mass distribution, tidal forces and torques ensue, and these can indeed cause chaotic motion, see, e.g. page 16 and 17 of these lecture notes.

In fact, I compared the trajectories of pointMass1.r_0[1] computed by Dymola and OpenModelica with Tolerance = 1e-6 and they are completely different after 7.5 seconds. Of course they start to diverge slightly (but enough to fail the CSV compare test) much earlier than that, exactly as in the case of Chua's circuit.

I then tried with Tolerance = 1e-8, and now they start to diverge slightly around time = 7.6 and dramatically around time = 12. If I ever saw a system showing tell-tale signs of chaotic motion, this is definitely one.

What is really interesting is also that the trajectories generated by Dymola with Tolerance 1e-6 and 1e-8 are completely different after time = 7.6. In fact, as @hubertus65 suggested to try out in some other ticket, even the trajectories generated by Dymola with the same Tolerance = 1e-6, using two different integration algorithms (Dassl and Radau II), are completely different after 7.5 s.

This proves quite compellingly that the motion is chaotic, i.e., any small perturbation is inevitably causing large arbitrarily large deviations if enough time elapses.

As to BevelGear, I have no rigorous proof of chaoticness at hand. However, the superposition of nonlinear oscillatory phenomena is a well-known path to chaotic motion, as in the case of a double-pendulum, or even as in the case of a single pendulum driven by a sinusoidal torque.

If I compare the trajectories of revolute1.axis.phi between Dymola and OpenModelica with Tolerance = 1e-6, they start diverging seriously around time = 0.5, while if I set Tolerance = 1e-10 they match very well until time = 0.72, then they diverge abruptly. Same behaviour when comparing simulations in Dymola with different tolerances and/or integration algorithms. This behaviour is also 100% consistent with chaotic motion.

We should always provide reference cases where the global numerical error with Tolerance = 1e-6 can be expected to be well below the tolerance of the CSV compare tool. For asymptotically stable trajectories, this is no big deal, because the effects of numerical errors tend to decay over time. For simply stable trajectories (e.g. an undamped oscillator driven by a sinusoidal force at the resonance frequency) you have to be a bit careful, because they tend to pile up over time. For chaotic trajectories, you need to stand clear from the edge of chaoticness with a robust margin, because the effects of arbitrarily small perturbations grow exponentially, until trajectories diverge dramatically from each other .

I think I proved beyond reasonable doubt that this the case for the three referenced MSL examples.

Of course one could question the StopTime's I am proposing, arguing that they are still too large. They are purely empirical and based on the criterion that the trajectories produced by Dymola and OpenModelica with Tolerance = 1e-6 should be closer than the CSV compare algorithm requires, plus some safety margin. If other tool vendors want them further reduced (within reason) I'm fine with that.

[maltelenz]

I am also in favor of shortening the reference results for models like these as proposed.

[casella]

@beutlich, I can very easily edit the CSV files and just remove the last lines past the required stopTime, then open a PR. This doesn't even require to re-run the simulations. Would that be fine for you?

[beutlich]

Well, my idea would be to rerun simulations for all models again (also where we removed/added comparison signals recently) and make it available with MSL 3.2.3+build.5.

[casella]

Sounds good.

In the future, I believe we should move towards a Continuous Integration approach. I hope to be able to work on that starting in September.

[AHaumer]

Another idea (different from having very short end time): change default parameters to obtain periodic behavior, and mention in the documentation which parameter values lead to chaotic behavior. See: #4477