Options of WhenEvent

Question

I'm using NDSolve and WhenEvent to simulate building heating systems and their discontinuous operation. Here is a simplified example.

The idea is quite simple. I have the following WhenEvents:

{WhenEvent[True && (tt[tau] <= 19) && (ta[tau] < 10), Q[tau] -> 10000],
WhenEvent[{tt[tau] >= 21, ta[tau] >= 10}, Q[tau] -> 0]}

I understand that there are some problems with WhenEvents if I use And as Head in the conditions. I found out that adding a True improves its behaviour.

The actual problem requires a large equation system. I don't want any event missed, because it causes singularity problems. According to the documentation, I should use

"IntegrateEvent" -> True, "DetectionMethod" -> "Interpolation"

Applying these into the differential equations

Cwirk = 50 25 3;
ta[tau_] := 10 - 10 Cos[Pi/(12) tau]
sol = Flatten[
       NDSolve[{Q[tau] - 200 (tt[tau] - ta[tau]) == Cwirk tt'[tau],

         WhenEvent[True && (tt[tau] <= 19) && (ta[tau] < 10), 
              Q[tau] -> 10000, "IntegrateEvent" -> True, 
              "DetectionMethod" -> "Interpolation"],

         WhenEvent[{tt[tau] >= 21, ta[tau] >= 10}, Q[tau] -> 0, 
          "IntegrateEvent" -> True, "DetectionMethod" -> "Interpolation"],

          tt[0] == 10, Q[0] == 10000}, {tt, Q}, {tau, 0, 24 7}, 
        DiscreteVariables -> Q]]

It gives the following errors

NDSolve::nlnum1: The function value {-11.,True} is not a list of numbers with dimensions {2} when the arguments are {0.,10.,10000.}.
Join::heads: Heads List and Experimental`NumericalFunction[{Hold[Function[{tau,tt,Q},{-21+tt,tt<=19&&Plus[<<2>>]<10}]],Apply},{0,{<<1>>}},<<2>>,{<<1>>},{Function[{NDSolve`Monitor$3956},Function[{NDSolve`Monitor$3956$3957,NDSolve`Monitor$3956$3958,NDSolve`Monitor$3956$3959},Internal`InheritedBlock[{Q,tau,tt},<<1>>]]],None,None}] at positions 1 and 2 are expected to be the same.
NDSolve::nlnum1: The function value {-11.,True} is not a list of numbers with dimensions {2} when the arguments are {0.,10.,10000.}.

Any ideas?

Answer 1

Events are processed in special ways

Consider the general form:

WhenEvent[ A && B , code; actions, options...]
          ^event ^

Though the title suggests options are a key interest in the OP, they are not the main problem underlying the error. One may put code before the actions, which are rules or keyword strings. I use it below, but it, too, is irrelevant to the problem.

I wish to focus on the event argument A && B. When the event has the head And, it is a special case that is treated in the following way, which leads many to misunderstand WhenEvent:

As an event in WhenEvent[], A && B is not logically equivalent to "A AND B"; rather, it is equivalent to "A IF B."

There are several special forms outlined in the docs for the principal part of the event A that get translated into appropriate event conditions. Expressions A that do not match the forms are treated as is and their T/F values are tracked. The gist of the docs is that an event is triggered when the value of A changes from False to True, provided the predicate B is True at the time.

An example of a special form is when A has the head Equal, such as f == 0. The equation f == 0 is virtually never True in a discrete, numerical solver, unless by some chance a numerical step lands f exactly equal to zero. So instead it is treated as equivalent to two events f > 0 and f < 0. A similar case is g == c, which is translated to the f == 0 case via f = g - c. The form Mod[..] == 0 is another special case. We'll come back to some of the special forms below, ones that seem related to the problem at hand.

A workaround for some cases

When "A AND B" is desired, one way is to hide the And in another function:

event[t_] := A && B; (* presumably A, B are expressions in terms of t *)
NDSolve[... WhenEvent[event[t], ...]...]

Now the head is event and not And. Further, event[t] is not evaluated symbolically to see if it results in a And[..] expression. In this case an event will be triggered when event[t] changes from False to True. This method can be used for other complicated logical combinations of conditions. (Beware: event[t] should not contain equalities, because they are virtually never True as discussed above.)

Alternatives or not:

One might try to come up with alternatives, especially you like to avoid defining unnecessary functions like event[]. Sometimes the obvious works and sometimes not:

If[A && B, True, False]  (* works *)
(A && B) == True         (* fails *)

I thought (A && B) == True was a clever way to hide the head And, but as noted above, the head Equal is another special case. This condition is transformed to

(A && B) - True == 0

(via the f = g - c transformation mentioned above). Internally a condition f == 0 is represented by a function that must have numerical values, which are used to detect a zero crossing. Obviously (A && B) - True is not numerical, and you get a NDSolve::nrnum1 error like the OP.

Example of various forms

To see how A && B compares with event[t] in all logical combinations, I'll use for A && B the following condition:

Sin[3 t] > 0 && Cos[t] > 0

The following will increment the discrete variable p at each A && B event, the variable q at each event[t] event, and as well, the variable r whenever A or B separately change from False to True.

Clear[event, p, q, r, x, t];
event[t_] := Sin[3 t] > 0 && Cos[t] > 0;
ss = NDSolve[{p[1] == 0, q[1] == 0, r[1] == 0,
    WhenEvent[Sin[3 t] > 0 && Cos[t] > 0, p[t] -> p[t] + 0.1],
    WhenEvent[event[t], q[t] -> q[t] + 0.1],
    WhenEvent[{Sin[3 t] > 0, Cos[t] > 0}, r[t] -> r[t] + 0.1],
    x'[t] == 1, x[1] == 1},
   {p, q, r}, {t, 1, 20},
   DiscreteVariables -> {p, q, r}, MaxStepSize -> 0.5];
Show[
 RegionPlot[{Sin[3 t] > 0, Cos[t] > 0}, {t, 1, 20}, {y, -1, 1.2}, 
  BoundaryStyle -> Thin], 
 Plot[{Sin[3 t], Cos[t], p[t] /. ss, q[t] /. ss, r[t] /. ss},
  {t, 1, 20}, PlotLegends -> {None, None, p, q, r}],
 AspectRatio -> 0.6
 ]

The shading from RegionPlot shows when A is true (blue) and B is true (golden).

A solution for the OP

Cwirk = 50 25 3;
ta[tau_] := 10 - 10 Cos[Pi/(12) tau];
evts = {};  (* optional: linked list for recording events; see plot *)
sol = Flatten[
   NDSolve[{Q[tau] - 200 (tt[tau] - ta[tau]) == Cwirk tt'[tau], 
     WhenEvent[{tt[tau] <= 19 && ta[tau] < 10, 
       ta[tau] < 10 && tt[tau] <= 19},
      evts = {evts, on[tau]}; Q[tau] -> 10000,
      "IntegrateEvent" -> True, "DetectionMethod" -> "Interpolation"],
      WhenEvent[{tt[tau] >= 21, ta[tau] >= 10},
      evts = {evts, off[tau]}; Q[tau] -> 0,
      "IntegrateEvent" -> True, "DetectionMethod" -> "Interpolation"],
      tt[0] == 10, Q[0] == 10000}, {tt, Q}, {tau, 0, 24 7}, 
    DiscreteVariables -> Q]];
evts = Flatten@evts;
Show[
 Plot[{tt[tau] /. sol, Q[tau]/1000 /. sol, ta[tau]}, {tau, 0, 24 7}, 
  PlotStyle -> {Automatic, Purple, Gray}, PlotPoints -> 500, 
  PlotLegends -> {tt, Q, ta}],
 Prolog -> {AbsoluteThickness[0.6],
   Red, InfiniteLine[{#, 0}, {0, 1}] & /@ Cases[evts, on[t_] :> t],
   Black, InfiniteLine[{#, 0}, {0, 1}] & /@ Cases[evts, off[t_] :> t]},
 PlotRange -> All, Frame -> True, AspectRatio -> 1/3, 
 AxesOrigin -> {0, 10}
 ]

(The red grid lines are events when the heat was turned on, and the black are when the heat was turned off.)

On the OP's problems and WhenEvent options

I have a guess, which could be interpreted as supported by the documentation, as to why the options lead to error message. From the tutorial Events and Discontinuities in Differential Equations, the description of the "DetectionMethod" option value is:

"Interpolation" -- interpolate the event function with a polynomial with the same order as the method.

I interpret (or guess) this means that the event needs to be defined by a numeric function derived from the condition A. When A is of the form f == 0, f > 0, etc., the expression f is used. When A && B is of the form True && ..., then A is True and this expression is used. Unfortunately True is not a numeric function and "DetectionMethod" -> "Interpolation" fails to work. Neither the error message nor the documentation clearly says this is not allowed, although I'm suggesting it can be inferred from the snippet above. In any case, it seems to me to be the case, which I will summarize in the following rule:

The option "DetectionMethod" -> "Interpolation" may be used in WhenEvent only when the event is defined by a numeric function (such as f == 0, f > 0, or f < 0).

(Side note: In the tutorial, it indicates "DetectionMethod" -> "Interpolation" automatically invokes "IntegrateEvent" -> True. So it is unnecessary to use both as above and in the OP.)