A couple of years ago I asked a question about solving n-body systems with NDSolve and detecting collisions with WhenEvent. I received an excellent answer, but unfortunately that code does not work anymore in Mathematica 10.2.0.

In particular, in v9 the integration continued after the collision (which occurs at t=2.39343) and properly updated the modified masses and velocities of particles, as can be seen in the resulting trajectories plot:

ParametricPlot3D[Evaluate[#@t & /@ pos], {t, 2.25, 2.43}, 
 PlotRange -> All, ImageSize -> 500, PlotPoints -> 500, 
 MaxRecursion -> 15, Boxed -> False, Axes -> False]

parametric plot from mma 9

In Mathematica 10, however, after the collision, the system falls apart and NDSolve can't produce sensible trajectories anymore - the parametric plot now looks like this:

enter image description here

How to fix the code in the linked answer for Mathematica 10 so that it will properly detect the collisions and update the masses and velocities of particles?

  • $\begingroup$ The first time I run the code, I get an error. It seems pids[[1]] is being evaluated before pids is defined. If I run the code a second time, it goes through without an error. The plot I get is the same as you V9 plot (V10.2, Mac OS): i.stack.imgur.com/J0oIL.png $\endgroup$ – Michael E2 Nov 24 '15 at 14:44
  • $\begingroup$ @MichaelE2 Interesting, I have the same specs. It happens on a fresh kernel. I get the errors too, but it doesn't affect the output. $\endgroup$ – shrx Nov 24 '15 at 14:56
  • 1
    $\begingroup$ No difference on a fresh kernel for me. If I put a random initialization, say, pids = {1, 2};, before whenevnt, it works fine. If I wrap the second argument of WhenEvent in Module[{pids},...], the error refers to pids$[[1]], which shows that NDSolve is inspecting the held code, since pids$ does not have a module number at the end. The error then gives a stiffness/singularity warning and integration stops at the collision. In that case, I get your V10 graph. Try the initialization dodge. See if it works. $\endgroup$ – Michael E2 Nov 24 '15 at 15:05
  • 2
    $\begingroup$ My guess is that when the part error on pids[[1]] happens, NDSolve decides to ignore the event. Seems like a bug to me. At the least, there ought to be a message warning the user that there was an error in the event and that the event is being ignored (assuming it is ignored and not merely misfiring). $\endgroup$ – Michael E2 Nov 24 '15 at 15:11
  • $\begingroup$ The error occurs in the NDSolve`ProcessEquations stage, which is not surprising. $\endgroup$ – Michael E2 Nov 24 '15 at 15:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.