Should Modelica.Constant.inf really be the largest representable number? #4479
Comments
I think the involvment of (practically unknown when modelling) mathematical manipulations will make it challenging to come up with a number for what "comfortably far" is, if that margin should ensure that the model does not run into floating-point arithmetic-related problems. Consider that one concrete example of problems discussed was the (numerically problematic) expression julia> 1e-13 + 1.7e308 + 1e-13 - (1e-13 + 1.7e308)
0.0Even if we reduce the large number to a value more "comfortably" away from DBL_MAX, this fails, because the floating-point number spacing at the large value is larger than the small number being added: julia> 1e-13 + 1e60 + 1e-13 - (1e-13 + 1e60)
0.0
julia> eps(1e60)
1.78405961588245e44So, in practice, I'm not sure what process to use (or even if there is one) to determine a number always "safe" from IEEE-754 problems. 🤷 Side remark on:
This is not accurate, the max representable does not have all julia> bitstring(1.7976931348623157E+308)
"0111111111101111111111111111111111111111111111111111111111111111" |
Agreed.
Or alternatively - every number is safe from IEEE-754 problems as the standard is well-designed. I would just close this as not an issue.
Obviously the sign-bit (first bit above I guess) isn't set for a positive number. However, I don't understand the position of the next zero, as it seems to be in the mantissa and that cannot be right. The max "exponent" (2047 for double) is reserved for actual Inf and NaN-variants in IEEE-754. Added: That bit-string is 75(?) characters long - double is 64-bits, and Intel's extended double is 80 bits. |
|
The more obvious question is: why is Otherwise, I recommend to rename In eFMI, we had to model all of this target independent, coming up with the following builtin functions:
That said, and as you see in above definitions, it is of no use to look into any bit patterns or considerations like programming language A does that and B this, since the truth comes from the hardware, where we can require it complies to IEEE 754. |
I must have made some mistake when copy/pasting from the terminal.... that bit string is indeed 64 bits long (and is correct in the terminal history): julia> typemax(Float64)
Inf
julia> prevfloat(typemax(Float64))
1.7976931348623157e308
julia> bitstring(prevfloat(typemax(Float64))) # to make sure I really got dbl_max
"0111111111101111111111111111111111111111111111111111111111111111"
julia> bitstring(1.7976931348623157e308)
"0111111111101111111111111111111111111111111111111111111111111111"and I see the same pattern in an online tool. Sorry about the confusion 😞 |
I agree in principle, but in practice what value should users use when the want a really large number that is numerically safe? That is a legitimate question that needs to have a proper answer. |
@bilderbuchi the smaller Constants.inf is not meant to make 1e-13 + inf - inf = 1e-13 reliable. It is meant to make sure that the division by zero is actually never computed, because it is in a branch that is only activated when time > very_very_very_large_value The problem is, the logic that handles time > very_very_very_large_value in OMC fails for some reason above 0.8e308, i.e. half of the largest real ever. Shall we spend time to fix it so that we can actually use 1.7976931348623157E+308 as very_very_very_large value? I don't think so. Who cares? There is no quantity in physics using SI units that ever goes that far. The diameter of the universe is 8.8e26 m, the age of the universe is 44e16 s, and the number of atoms in the universe is 1e80. Why bother getting all the way to 1.7976931348623157E+308? I don't get it. Of course an easy solution for us is to change the value of inf in ModelicaServices (which is tool dependent) and that's the end of the story. But I think it would be best to define a "safe" very_very_very_large value in Modelica.Constants. My take is that the previous 1e60 was not too bad: it has the order of magnitude of the fifth root of 1.7e308, and still largely exceed any SI unit value we may ever find in a Modelica model. As I wrote, we should probably not call it In any case, referring to #4478, if we keep Modelica.Constant.inf at 1.7e308 (which is not numerically safe by definition, since it is a limit value) , I guess we should then put something more reasonable as the default value for V_pulse.PW and V_pulse.PER. That model is meant for circuit simulation, so there is absolutely no reason to get to such large values to mean "infinity". |
Yes, OpenModelica organization should fix that. |
|
Just in case one wants to run the circuit simulation for 1.2e308 seconds... Seriously, the reason why it shouldn't is explained in #4478 |
Because AFAIK +INF is defined by IEEE 754 but not by the Modelica Language Specification. Ditto for NaN |
As an aside, maybe that should be amended, too? I feel like more compatibility with IEEE-754 on a language level can only be a good thing, considering the domain Modelica operates in. |
|
Does anybody at OSMC understand why OpenModelica has a problem for values larger than approximately inf/2? |
The problem is that equation-reasoning and non-finite numbers clash. The minor issue is that we generate an error for 1/0 or sqrt(-2), but with IEEE-754 they would work - but either we still generate an error if used for the derivative (etc), or we generate a simulation result that is mostly NaN. I'm not saying that it would be impossible to have them (especially in some restricted form), just that it would require a lot more work. One possibility would be to have some other data-type for it, similarly as Complex numbers allow Modelica.ComplexMath.sqrt(-2). |
Yes, see my comment in #4478. In a nutshell, you can't reliably add a positive offset to a number which is already the maximum possible representable value. It is as simple as that. It's not a problem with OpenModelica. |
That is contradicted by your previous statement that it works to about 0.8e308 and fails for 0.9e308 (and for 1.7...e308). I'm sorry to say, but this is clearly an OpenModelica problem not a MAP-Lib issue. If you want to have a patched version of ModelicaServices to work around that in Openmodelica it is one thing - but it is still not a MAP-Lib issue. |
@bilderbuchi, motion of order. This issue tracker is about the Modelica Standard Library, we shouldn't start discussing a-sides here, otherwise managing these issues becomes very diffcult. Changing the semantics of the Modelica language by adding Inf and NaN is a big deal, which could be discussed in the proper place, e.g. a ticket in https://github.com/modelica/ModelicaSpecification/issues. As I see, there is absolutely no need of doing something this disruptive to solve the issue of this ticket, which is, I recall, about a regression of MSL that should be resolved in a matter of a few weeks. |
@HansOlsson I can't see why it should. My statment is that adding any positive number to Anyway, I think I can close this ticket, after all these discussions I made up my mind and I am happy at keeping I am sorry if it took a while to figure this out, I guess everybody understands that this was not a trivial problem at all |
If OpenModelica fails when adding a small value to 0.9e308 (as you stated) I don't see that we can use the fact that OpenModelica cannot handle adding a small value to 1.7...e308 to confirm that adding a small value to 1.7...e308 is an unreliable operation in general. |
Absolutely. My argument is exactly the opposite: adding any value (small or large) to I would also add that I don't see much value in discussing forever what happens when we make operations with numbers of the order of magnitude of 1e308, in a declarative modelling environment. I would find a discussion about how many angels can dance on the pin of a needle more useful. What's the point of doing that? Ultimately, who cares? We should avoid those operations in the first place in our models. |
|
I was asked to speak my verdict on this matter. My judgement on the particular matter:
Verdict: The model with its default parameterization is invalid as soon as T0 > 0. Underlying root cause The underlying problem is the following:
The combination of these two points is however simply asking for trouble. It is either or. If you provide a maximum floating point number, you somehow have to define the behavior around that number because people will start (ab-)using it. The models above are a wonderful example for this. If you choose not to define the behaviour, you shall not provide this number in the first place. Suggested short term remmedy Since this problem only seems to occur in OpenModelica, I think we can turn a blind eye on this and ask OM to set Suggested permanent remmedy Since we provide
A defintiion in this style would still leave some freedom for the actual numerical implementation. However this is just a quick suggestion. Not entirely thought through... I would neet to check whether this actually the case for IEEE 754 We should also rename it to Final comment Technically, my suggested solution is incorrect, since MSL 4.1 shall be based on Modelica 3.6 and not on 3.7 or whatever. However, I have the feeling that the usage of Edit |
|
Thought a bit more on this matter.... I think that there are good reasons not to define the floating point behavior in the Modelica language If we do not define the behavior, then the specifics of a floating point arithmetic have no place in Modelica.Constants We have to define such constants then by intent and way of usage. For In absence of true infinity,
This would clarify, what computations can be expected to run. |
In fact, it is also invalid for |
|
@dzimmer: You following assessment is off, I am afraid:
The underlying, general issue is that any saturated maths scheme violates basic laws of algebra, like commutativity, distributivity etc. As a simple example, assume I also don't buy the argument that Modelica should not just comply to IEEE 754. What kind of magic hardware are you guys using? Did anybody in this thread ever encounter an application of Modelica simulation not on IEEE 754 compliant hardware? I can imagine such case for some exotic embedded systems, but that is not the current hardware people run Modelica IDEs or simulations on. Finally, I also have the feeling that |
Whether you buy it or not, there is no mention of IEEE 754 in the language spec, except for the procedural external functions, because the language is declarative. Of course the code generated by a Modelica compiler in the 21st century eventually complies to IEEE 754. But that doesn't mean the original Modelica code does, nor that it should. Modelica in general is about writing equations, not computer code. Consider the following MWE package MWE
model M
Modelica.Blocks.Interfaces.RealInput u(nominal = 1e6);
Real y;
parameter Real p = Modelica.Constants.inf;
equation
y = p - u;
end M;
model S
M m;
parameter Real r(min = Modelica.Constants.small);
equation
r*m.y = time;
end S;
end MWE;The symbolic solution of the equations of which according to IEEE 754 should be OK without underflow or overflow. However, there is no guarantee that the solution is actually computed in that way and, in general, that a symbolic solution can be computed at all, or that it is the most efficient way to get a solution. A tool could (at least in principle), exploit the nominal attribute of I could come up with more contrived MWEs, but to me it is clear that in general one cannot make statements about the IEEE 754 validity of individual Modelica models by looking at their code, because of the a-causal nature of the language. It may still be possible to introduce NaN, Inf and the associated semantics in Modelica, but making sure that it is preserved through the code generation process will definitely be a non-trivial task and may lead to significant overhead. It could be an interesting research question to establish whether this is doable and worth the (compile-time and run-time) cost.
I totally agree with that, and in fact the solution that I proposed in #4478 avoids using |
|
I agree that if one chooses to specify the floating point behavior then it is best to simpl comply to IEEE 754. My proposal was not thought out. Forget about that. However, the current Modelica Language does not specify IEEE 754. It even explicitly defines the behavior in case of overflows to be undefined. (Sec 3.3 $2) What have been historic reasons for that I do not know. Francesco and me speculated that is because equation-based models are declarative in nature and may be processed differently. My problem is that when you do not specify the floating point behaviour then you should not provide constants that are specific (edge) cases of floating point numerics. This invevitable will lead to problems.
I am basically fine with both options. The current combination is just unlucky... |
No, it is not at all that bad. For In order to run into the undefined behavior of an overflow, one needs to add pretty big numbers to |
This would be OK assuming that the expressions in the Modelica code are are always interpreted literally. The problem is, they aren't. A Modelica tool is free to apply transformations to expressions and equations that are equivalent from a mathematical point of view, but not from a numerical point of view. For example, it can rescale variables and equations based on nominal attributes, which is in general a good idea if SI units are used, see DOI:10.1145/3158191.3158192, in order to obtain better scaled equations for numerical solvers. Or replace expressions with other expressions that are mathematically equivalent but, hopefully, better numerically behaved. Or solve/simplify equations and inequalities by symbolic manipulation. Considering the case that triggered this discussion, see #4478, the simulation has I'm not saying that doing this time scaling is a good strategy for simulation nor that tools should actually implement it, this is just an example. What I am saying is that we should take the fact that Modelica is a declarative mathematical modelling language seriously when discussing this matter. |
Follow up of #4478.
The discussion in #2056 ended up in PR #4042 for MSL 4.1.0, in particular:
ModelicaStandardLibrary/ModelicaServices/package.mo
Lines 184 to 185 in 8e7876b
This is the actual largest number that can be represented with a double-precision number, i.e. a number with all ones in all bits. Change it up by a small fraction and it will overflow. Similar arguments can be taken for eps.
Now, regardless what we write in the comments we put on those definitions,
Modelica.Constants.infis used by library developers (including MSL) with the meaning of "very large number, approximating infinity from all practical purposes". That's a fact of life that we need to live with, unless we give them other options. BTW, a library developer will use such a constant with more confidence than a made up large value, since it is defined in the MSL.For example, if I have a system which triggers an event at time
T, I can setT = Modelica.Constants.infas a default, meaning "never" in practice. This means you have code such asso we are using
Modelica.Constants.infas a regular numerical value for end usage in equation-based Modelica models, not for some low-level algorithmic programming task.Now, if Modelica was an imperative language, this would be OK, because such a value would be interpreted literally. But we know very well that Modelica is no such thing: with the exception of algorithms in functions, every expression can be subject to mathematical manipulation: scaling with nominal attributes, symbolic solution, symbolic differentiation, zero-crossing function generation, index reduction, etc. etc.), including those that contain
Modelica.Constants.inf. What if the code for event generation requires to handle values slightly larger than the RHS of the inequality, for any reason? So, the question is, what is the "largest representable number" for a declarative modelling language as Modelica? Is it really the borderline value with all bits set to one?This is not a purely hypothetical scenario, but in fact happened in OpenModelica after #4042 was merged in, see OpenModelica/OpenModelica#13010: Spice3 V_pulse generators have default values of period and pulse length both set to Modelica.Constants.inf. It turns out that if values up to half of the largest representable double-precision point are used for that constant, everything works fine, while if the extreme value is used, something breaks with the event management process, most likely because of overflow issues.
Maybe a good solution to this issue could be to keep
Modelica.Constants.infdefined as "the largest possible value" (just in case it is needed inside some complicated algorithm), but then define something likeModelica.Constants.largethat can be considered as infinite for all practical modelling purposes, but is still comfortably far from the numerical overflow limit. This should be used in all cases where the constant is not meant to be used literally in an algorithm, but possibly subject to some kind of symbolic/numeric manipulationWe should then change all uses of Modelica.Constants.inf in the MSL that are just meant to use a very large value to Modelica.Constants.large.
Alternatively, I would seriously consider setting Modelica.Constants.inf back to 1e60.
I'd like to get some feedback also from library developers on this topic. @hubertus65, @GallLeo, @mwetter, @wischhusen, @AlesVojacek, @bilderbuchi, @adriguir, feel free to comment.
The text was updated successfully, but these errors were encountered: