# Weird behaviour in WhenEvent

Consider the following code (taken from question 158857):

g = 9.81;
eqs = {x1''[t] + g == 0, x2''[t] + g == 0};
ci = {x1[0] == 1, x2[0] == 2, x1'[0] == 0, x2'[0] == 0};
events := {WhenEvent[x1[t] == 0, x1'[t] -> -x1'[t]],
WhenEvent[x2[t] == x1[t], v1 = x1'[t]; x1'[t] -> x2'[t]],
WhenEvent[x2[t] == x1[t], (* -> *) 1; (* <- *) x2'[t] -> v1]}
sol = NDSolveValue[eqs~Join~ci~Join~events, {x1, x2}, {t, 0, 5}]

Plot[{sol[[1]][t], sol[[2]][t]}, {t, 0, 5}]


It returns

The first WhenEvent accounts for the bouncing on the ground ($y=0$) while the other two aims at exchanging the derivatives.

Now if I remove the 1; in the last WhenEvent, I get the following error

NDSolveValue::nlnum1: The function value {v1} is not a list of numbers with dimensions {1} when the arguments are {0.790166,0.9375,1.10736,0.9375,-7.75153}.

Why is only working if I add something before x1'[t] -> x2'[t] ? I tried on MMA 11.0 and 10.3.

• Module[{}, x2'[t] -> v1] makes it work as well. Odd. Commented Oct 31, 2017 at 18:25
• you can also do this more cleanly with a single WhenEvent : {x1'[t], x2'[t]} -> {x2'[t], x1'[t]}] ( not really answering the question though ) Commented Oct 31, 2017 at 18:28
• @george2079 1. Strange. 2. I thought that was the first thing I had tried. So much simpler... Thank you! Commented Oct 31, 2017 at 18:32
• I suspect that the order that the two identical events are dealt with differs. If I insert a Print statement into the first x2[t]==x1[t] event, that is printed after the error message (so in fact, v1 hasn't been defined yet). Commented Oct 31, 2017 at 21:05
• In v9.0.1 there's no warning, but the result is wrong: i.sstatic.net/lUQGO.png Commented Nov 1, 2017 at 3:08

I hope my answer does not come too late. You should include v1 into DiscreteVariables

g = 9.81;
eqs = {x1''[t] + g == 0, x2''[t] + g == 0};
ci = {x1[0] == 1, x2[0] == 2, x1'[0] == 0, x2'[0] == 0, v1[0] == 0};
events := {WhenEvent[x1[t] == 0, x1'[t] -> -x1'[t]],
WhenEvent[
x2[t] == x1[t], {v1[t] -> x1'[t], x1'[t] -> x2'[t],
x2'[t] -> v1[t]}]}
sol = NDSolveValue[eqs~Join~ci~Join~events, {x1, x2, v1}, {t, 0, 5},
DiscreteVariables -> v1]
Plot[{sol[[1]][t], sol[[2]][t]}, {t, 0, 5}, PlotRange -> All]