Dear Colleagues,
I need to invert a very complicated function (i.e. given y = f(x), I wish to find x = g(y)) and neither the inverse nor the function can be given analytically. I have written a simple test program (attached) that uses GetX() to invert an equation of the form y = exp(x). It appears to work in root version 5.34 but not version 6.18.
Can anyone kindly explain what I have to do to make it work in version 6.18? shine.c (1.1 KB)
Thanks @vpadulan for the very prompt reply. Is there any chance that you can compile the code on root v. 5.34? In that version there is no error message and it does find the inverse function. Also, by inspection, the equation ‘0 = exp(x)-y’ must contain the root ‘x = ln(y)’. The question then is why does root v. 6.18 not find it.
It seems to me that the there is no x in which f(x) = 0, that means that the errors you are getting in more recent versions of ROOT are correctly reporting this fact. I’ve modified the call to GetX in your root_int to GetX(1,0,10) and the programs runs without errors outputting the following figure
Thanks for looking at it again @vpadulan. I believe, however, that a significant point is being overlooked. If the original function is y = exp(x), then the inverse is x = ln(y) and for x = 0, y = 1. In other words the curve should go through the point ( x=0, y = 1) as it does using the original equation I ran in root v. 5.34.
Furthermore, if we change the function to another one, say, y = sin(x), the same problem still occurs with root v. 6.18, while the solution given by root 5.34 is correct.
Does anybody know what I can do to make the program work in root 6.18?
As far as I can tell it has something to do with the function. If I choose a polynomial such as y = pow(x,2), then root 6.18 has no problem inverting it. If, however, I choose a function that should be in the standard library e.g. y = sin(x), y = exp(x), y = cos(x) etc, then it fails.
Not an expert on this, but I had a quick look in the code. (Sorry, I don’t have time to check the math this week.)
It looks like the BrentRootFinder, which is used to solve the GetX problem doesn’t find a clear minimum, that is, there is either a degeneracy or as @vpadulan pointed out, there is simply no value where f(x) = 0 .
To test, I hacked the function that complains that there’s no inverse in ROOT 6.22. It will now just return zero whenever there’s no inverse value, and the result is this:
To be clear: At all points that are zero in this plot, the function may actually be undefined.
I see two alternatives:
In 6.18, ROOT could be trying to invert the wrong function. If you think that’s the case, I guess we have to wait for @moneta to return from vacation. Hopefully, my plot can show if the function it’s trying to invert is correct.
If it’s inverting the correct function, there’s an ambiguity or the function never reaches zero. In that case, the definition range of the original has to be altered such that it is invertible everywhere.
Part of what you say makes a lot of sense, @StephanH:that BrentRootFinder, for whatever reason, doesn’t find a clear minimum. It’s wierd though that it can find it in root v. 5.34, but not in root v. 6.18.
Of the two alternatives you offered, I ran a test program that I think eliminates (2). y = sin(x)reaches 0 and is not ambiguous on the range GetX(0, 0, 0.6). However, the problem still persists in root v. 6.18.
I was told that there was indeed a change in the BrentRootFinder:
In ROOT 5, it would return a degenerate minimum as a valid minimum. In ROOT 6, it will warn when the function cannot be inverted.
This indicates that in your case, the problem is indeed that the minimum is degenerate. Could you verify this?
I have checked and the problem is in your function. You try to find a root from exp(-x)-p. When p is less than 1 there is no root. What the code returns is correct.
5.34 was probably havinga bug and probably returning a wrong result, so forget about that version
Again it is a problem in the function.
For sin(x)-y the root is for x = arcsin(y) and it is correct, but the interval should be [0,1] and not [0, 0.7]
For example for y = 0.69 you have the root for x ~ 0.76 that is larger than 0.7.
Again the code is correct to signal that there is no root. You should catch that condition (GetX returns a NaN because there is no root) and act in your code
You were right @moneta. It seems though that Root 6 has less tolerance for the kind of loose coding I had in the original versions of my program. The answers from Root 5 are correct away from the singularities and at the singularities the code just returns default values.