0

I build a simple model to understand the concept of "Discrete expressions", here is the code:

model Trywhen
  parameter Real B[ :] =   {1.0, 2.0, 3.0};
algorithm 
    when time>=0.5 then
      Modelica.Utilities.Streams.print("message");
    end when;
  annotation (uses(Modelica(version="3.2.3")));
end Trywhen;

But when checking the model, I got an error showing that "time==0.5" isn't a discrete expression. enter image description here

enter image description here

If I change time==0.5 to time>=0.5, the model would pass the check.

enter image description here

And if I use if-clause to when-clause, the model works fine but with a warning showing that "Variables of type Real cannot be compared for equality."

enter image description here

My questions are:

  1. Why time==0.5 is NOT a discrete expression?
  2. Why Variables of type Real cannot be compared for equality? It seems common when comparing two variables of type Real.
Jack
  • 1,094
  • 6
  • 16

1 Answers1

5

The first question is not important, since time==0.5 is not allowed.

The second question is the important one: Comparing reals for equality is common in other languages, and also a common source of errors - unless special care is taken.

Merely using the floating point compare of the processor is a really bad idea on idea on some processors (like Intel) that mix 80-bit and 64-bit floating point numbers (or comes with a performance penalty), and also in other cases it may not work as intended. In this case 0.5 can be represented as a floating point number, but 0.1 and 0.2 cannot.

Often abs(x-y)<eps is a good alternative, but it depends on the intended use and the eps depends on additional factors; not only machine precision but also which algorithm is used to compute x and y and its error propagation.

In Modelica the problems are worse than in many other languages, since tools are allowed to optimize expressions a lot more (including symbolic manipulations) - which makes it even harder to figure out a good value for eps.

All those problems mean that it was decided to not allow comparison for equality - and require something more appropriate.

In particular if you know that you will only approach equality from one direction you can avoid many of the problems. In this case time is increasing, so if it has been >0.5 at an event it will not be <=0.5 at a later event, and when will only trigger the first time the expression becomes true.

Therefore when time>=0.5 will only trigger once, and will trigger about when time==0.5, so it is a good alternative. However, there might be some numerical inaccuracies and thus it might trigger at 0.500000000000001.

Hans Olsson
  • 11,123
  • 15
  • 38
  • But `x>=y` still needs to be transformed into the form of "x-y>=eps", so the computer could handle it, right? So how should the machine determine the value of eps when using `x>=y`? – Jack Sep 12 '20 at 09:20
  • If it is not optimized into a time event (where we know the exact time it triggers: 5.000), a numerical solver will search for the first time where the condition triggers; this solver will contain its own epsilon value and is probably related to how large time steps you take and what would be considered a tiny value at this point in time. In OpenModelica, the maximum number of iterations is currently `1 + ceil(log(fabs(*b - *a)/(MINIMAL_STEP_SIZE + MINIMAL_STEP_SIZE*fabs(*b-*a)))/log(2))` where b and a are the times to look for zero crossings between (unless the user overrides it)... – sjoelund.se Sep 14 '20 at 07:56
  • For x>=y we can view it as a binary choice where we want to stay in the previous state (for stability), and the choice of epsilon isn't that critical. One implication is that after changing state then x>=y is (almost certainly) actually correct . That part we cannot guarantee for x==y. – Hans Olsson Sep 14 '20 at 09:32