I Runge-Kutta 4 w/ some sugar on the top: How to do error approximation?

[bremenfallturm]

Decimals are written as X.XXX... and not X, XXX... It's just a nit, don't mean to make it a big deal.

Aha, I see what you mean. I copied it from where I have noted the data and it ended up like that. Thanks for pointing it out!

 

[BvU]

The assignment is asking for "safe digits", so a relative error in my case. Here is an example of how I've worked through the first question (the one that just involves solving the differential equation):
Remember it asks for the coordinates of the first two bounces!

The assignment text is confusing me (again). Your 'lots of numbers' and the $\sum|\Delta|$ propagation method give me 6% error in bounce 1, 7% in bounce 2 and 5 in $y_{x=0}$. Huge ! (and the $\sqrt {\sum\Delta^2}$ about 3.5%, also huge).
Note that the uncertainty in the error estimate is considerable and only seldom warrants giving more than 1.5 digit !
(one digit, two if the first is a 1) .

"The answers should have at least 6-8 correct digits, somewhat depending on question" (yes, this is what the assignment says)
You also need to assume that every given value has an uncertainty of 1%. Use that to approximate the tabular error, and then include that in the calculation of te total error of your program.

So as an answer I can say that the uncertainty of bounce 1 occurs at −0,3534969571 ± 6,08%.

Which comes down to 0.35 $\pm$ 0.02
Two digits, no more.

What's the point of doing a 6-8 correct digits calculation if these uncertainties in the parameters are so big ?

[edit] can you check your 'lots of numbers' ? I find it strange that the first bounce occurs further to the left if $k_x$ is lowered ....




$\$

 

[bremenfallturm]

What's the point of doing a 6-8 correct digits calculation if these uncertainties in the parameters are so big ?

I see what you mean, it has of course frustrated me a bit too. We're supposed to consider two scenarios: the first one assuming that the input data is perfect (this yields 6-8 safe digits) and the second one assuming that the input data has that $\pm 1\%$ error (this yields 2 safe digits as we both concluded).
I will just assume that they want the first scenario to have the requested accuracy and the second one is separate to show that I have studied how errors affect the assignment.
I think the high requirement for precision is to show that we know how to apply efficient methods and predict errors, though it's not clearly stated which scenario that should have the requested accuracy, only that a "total error estimate" should be included.

 

[pbuk]

I will just assume that they want the first scenario to have the requested accuracy and the second one is separate to show that I have studied how errors affect the assignment.

Yes, this is how I interpret it too.

 

[bremenfallturm]

Thank you guys! So now I'm actually finalizing the thing that the topic title asks about: I'm putting together an answer!
Once again, please forgive me for not knowing these things that probably seem obvious to you. (Reminding you again that I'm a first year student and haven't been properly taught about error approximation!)

So let's say I want to provide an answer to everything.

1. Some answers have more accurate digits (7-9), does it make sense to round all numerical solutions to 6 digits or should I include all of them? (is it good practice to answer with all values having the same number of digits or should I not bother).

2. I am asked for a relative error, does it make sense to indicate the error boundary in percentages like this: $$\pm 0.5\cdot 10^{-6}\%$$?

3. I am also asked to indicate the presentation error. That is the error caused by rouding (for example, the presentation error of $1,2345$ is $\underbrace{1,2345689}_{\text{non-rounded value}}-1,2345=0,0000689$ if I round to 5 digits. Should I take that value, call it $E_{presentation}$ and express it in percentage like $E_{presentation pct}=\frac{E_{presentation}}{\text{non-rounded value}}$, or should I not bother about expressing it in percentage?

 

[bremenfallturm]

I double-checked my table, looks like you're right. k_x down gives −0,352986712689222.
I have another question. I have from above an uncertainty that is >5%. Doesn't that mean that no digits are safe? Since 0,5e-1 implies once safe digit. What do I answer in that case?

 

[BvU]

I double-checked my table,

Did you ? And did you notice something strange with with $k_y$ too ?

looks like you're right. k_x down gives −0,352986712689222.

I get -0.352986712689222
and for $k_{x, up}$ I get -0.35521856 and not -0.352986712689222

I have another question. I have from above an uncertainty that is >5%. Doesn't that mean that no digits are safe? Since 0,5e-1 implies once safe digit. What do I answer in that case?

We've whittled it down to 4% I think... But even 6% warrants two digits: 0.353 $\pm$ 5% is 0.353 $\pm$ 0.018 (the 1.5 digit if you firmly believe the accuracy if the error) or at worst 0.35 $\pm$ 0.2

And I think adding absolute $\Delta$ is way too pessimistic; adding squares and taking the root is much more sensible.

Also see #52

will look at #55 tomorrow. But $\pm 0.5\cdot 10^{-6}\%$ must be a mistake

$\$

 

[bremenfallturm]

Did you ? And did you notice something strange with with $k_y$ too ?


I get -0.352986712689222
and for $k_{x, up}$ I get -0.35521856 and not -0.352986712689222



We've whittled it down to 4% I think... But even 6% warrants two digits: 0.353 $\pm$ 5% is 0.353 $\pm$ 0.018 (the 1.5 digit if you firmly believe the accuracy if the error) or at worst 0.35 $\pm$ 0.2

And I think adding absolute $\Delta$ is way too pessimistic; adding squares and taking the root is much more sensible.

Also see #52

will look at #55 tomorrow. But $\pm 0.5\cdot 10^{-6}\%$ must be a mistake

$\$

Aha, I see you. Embarrasing, sorry!
I've always thought that a relative error of $0.5\cdot10^{-n}$ where $n$ is the number implies $n$ safe digits. Is that wrong?
It looks like it's the change in velocity that is causing such huge error.

 

[pbuk]

I've always thought that a relative error of $0.5\cdot10^{-n}$ where $n$ is the number implies $n$ safe digits. Is that wrong?

Yes it is wrong - that is an absolute error.

Edit: it could also be a relative error - see #61.

 

[pbuk]

And I think adding absolute $\Delta$ is way too pessimistic; adding squares and taking the root is much more sensible.

This depends on context - I would be inclined to do both, stating the assumption that must be made in order to use the second approach.

 

[bremenfallturm]

Yes it is wrong - that is an absolute error.

Aha, thank you! So a relative error of $0.5\cdot 10^{-n}$ implies $n+1$ safe digits? According to
https://math.stackexchange.com/ques...-relative-error-give-number-of-correct-digits . (that's good - or at least better - news!)

This depends on context - I would be inclined to do both, stating the assumption that must be made in order to use the second approach.

My lecturer uses $\sum \left|\Delta \right|$ when using the method I described, i.e. not the propagation suggested by BvU. I assume that the course wants me to use that one, even if it is pessimistic.

 

[pbuk]

Aha, thank you! So a relative error of $0.5\cdot 10^{-n}$ implies $n+1$ safe digits? According to
https://math.stackexchange.com/ques...-relative-error-give-number-of-correct-digits . (that's good - or at least better - news!)

Hmmm, there is the danger of confusion here. Let's make sure we are all clear on the definitions:

If the true value of a quantity is $x$ and the measured value is $x_0$ then the absolute error is $x_0 - x$ and the relative error is $\frac x{x_0} -1$.​

https://mathworld.wolfram.com/AbsoluteError.html
https://mathworld.wolfram.com/RelativeError.html

So $0.5\times10^{-n}$ could be either a relative error or an absolute error.

My lecturer uses $\sum \left|\Delta \right|$ when using the method I described, i.e. not the propagation suggested by BvU. I assume that the course wants me to use that one, even if it is pessimistic.

I agree with your lecturer, and it is not pessimistic: in the absence of other information it is the right method to use. The method suggested by BvU is optimistic because it assumes that the errors are independent random variables and in this problem there is no reason to indicate that that is true - which is why I suggested that if you did use it you should state that assumption.

 

[bremenfallturm]

Hello, sorry for not responding here sooner. I have waited for this project to be reviewed. Unfortunately, I got some adjustments to make - not because of the errors (so I guess we're all good on our definitions).

But unfortunately, using the current algorithm to detect zero crossing is "too inefficient, do interpolation instead" to quote the feedback I got. Here is the current algorithm:
Remember the problem is that I need to flip the sign of the y derivative at y=0.

1. Run Runge-Kutta with a "coarse" step length (0.5e-2) until a y solution crosses 0.
2. Use the secant (bisection) method to find a step length that will make us land at exactly y=0, with machine epsilon error acceptable.
3. Reset step length, change sign of derivative and continue iterating.

So, now I'm trying to figure out how to "do interpolation instead" to find the zero crossings.

It doesn't really make sense to me.

Let's say I find the 2 t (time) points closest to 0 and perform linear interpolation between them. I can then interpolate what time t that y is equal to 0 at, but the Runge-Kutta solutions vector is dependent on multiple things: $x, y, \dot x, \dot y$.
I fail to see what to do after finding this t. Do I reset to the previous solution, flip the derivative, and then recalculate with the interpolated t value? If so, I do not get a good precision at all and the graph looks really wonky:

What is meant by "do interpolation instead"?
Should I do interpolation over each of these: $x, y, \dot x, \dot y$?

 

[bremenfallturm]

Edit: I see what is meant, I interpolated over everything that the Runge-Kutta 4 solution vector is dependent on. Results look promising. I am tired and it's late, will elaborate further tomorrow.

 

[bremenfallturm]

Hello, time for an update I guess!
I hope I can get some last help from this topic, after that I'm hopefully done with the project.
I decided to do interpolation to find when y=0 the following way:

1. When y<0, solve for one more iteration so we have 2 "grid points" below y=0.
2. Take those 2 grid points and the 2 last grid points above y=0. Create an interpolation polynomial (of degree 3) between these points and interpolate $y, \dot y, x, \dot x$ as functions of $t$.
3. Solve using the bisection method for $y_{interpolated}(t)=0$.
4. Using the value obtained in (3), let's call it $t_0$, plug it into all the interpolation polynomials ($y_{interpolated}(t_0), \dot y_{interpolated}(t_0), x_{interpolated}(t_0), \dot x_{interpolated}(t_0)$. $\dot y$ changes sign according to the assignment.
5. The solutions obtained in (4) is set as being the last solution. All solutions below y=0 are scrapped, and Runge-Kutta solves for the next values using $t_0$ and the interpolated values as the last true solution

This does work, so I guess it's fine! However, I do have some new questions which also are (sigh) related to error calculation:
1. I am asked to have $6-8$ safe *digits*, aka. relative error which we have concluded in the past. With the interpolation around $y=0$ comes a new error, but it doesn't make sense to me how to calculate the relative error. When approximating the interpolation error, I compare my result with a fourth degree polynomial evalulated at the same point.
Here is an example for the two values of $y_{interpolated}(t_0)$ (see above for definition of $t_0$)
Degree 3: -6.938893903907228e-18
Degree 4: 2.081668171172169e-17
So the absolute error in this case is $\left| 6,938893903907228\cdot 10^{-18} - 2,081668171172169 \cdot 10^{-17} \right| \approx 2,8\cdot 10^{-17}$
But if I want to calculate the relative error, I take the absolute error divided by $6,938893903907228\cdot 10^{-18}$ which will turn out huge.
In some other cases, the value of the interpolation polynomial at degree 3 will turn out as $0$, which is good because I want to find when $y=0$. However, if I try to calculate the relative error in that case, I will get division by zero.
Since I use relative errors for everything else, I want to figure out how to calculate it for the interpolation polynomials.
2. Bonus question The interpolation around $y=0$ includes a small difference between the $t$ values when a fine step size is being used. For example, to solve the assignment which asks for a speed, to get the wanted accuracy, I have to interpolate between values that differ quite little in magnitude: here is an example:
2.685156250000135e-01 2.685937500000135e-01 2.686718750000135e-01 2.687500000000135e
The method does provide the same result as with the old algorithm, but I'm still a little worried regarding the interpolation accuracy - let's say I use these datapoints to create a coefficient matrix for the interpolation polynomial. Will that matrix not have quite a bad condition number?
I guess if it works it works, but is there any way I can improve the condition number (it has only been briefly mentioned in my course)

I hope this is clear and you can follow.

 

[pbuk]

In some other cases, the value of the interpolation polynomial at degree 3 will turn out as $0$, which is good because I want to find when $y=0$.

Yes: you want to find when $y = 0$, so you are looking for 8 digits of accuracy in $t$, not ## y #!

2. I have to interpolate between values that differ quite little in magnitude: here is an example:
2.685156250000135e-01 2.685937500000135e-01 2.686718750000135e-01 2.687500000000135e

Those numbers differ by a factor of approximately 1 + 2.9e-4; Machine epsilon is 1.1e-16 so you have about 16-4 = 12 digits of accuracy which is plenty.