Absolute and relative tolerance definitions

I am trying to understand the following Matlab definitions for RelTol and AbsTol parameters:

RelTol — This tolerance is a measure of the error relative to the size of each solution component. Roughly, it controls the number of correct digits in all solution components, except those smaller than thresholds AbsTol(i). 

The default, 1e-3, corresponds to 0.1% accuracy.

AbsTol — AbsTol(i) is a threshold below which the value of the ith solution component is unimportant. The absolute error tolerances determine the accuracy when the solution approaches zero.

I do not understand why AbsTol determines the accuracy when the solution approaches zero (indeed, if the solution of my problem is a circular orbit of 7000 km radius this does not meet) and why RelTol controls the number of correct digits in all solution components, except those smaller than thresholds AbsTol(i). I would like to have simpler and understandable definitions.

Answer 1

[Jan]

The relative tolerance between two values X and Y is:

abs(X - Y) / min(abs(X), abs(Y))

But if X or Y becomes (near to) zero, the tolerance grows until infinity. Then the absolute tolerance is a better limit:

abs(X - Y)

[EDITED] The relative error "abs(X - Y) / min(abs(X), abs(Y))" determines the number of "equal digits" between X and Y, if X and Y are >> 1. Just try it:

abs(7000 - 6000) / min(abs(7000), abs(6000)) < 0.1  % FALSE
abs(7000 - 6900) / min(abs(7000), abs(6900)) < 0.1  % TRUE
abs(7000 - 6900) / min(abs(7000), abs(6900)) < 0.01 % FALSE

The term "equal digits" is not meant very strict here. This is more a rule of thumb.

But if you get near to zero, this measurement fails, because the division by a small number let explode distances. Then the absolute distance is applied as limit, because it is more likely to match the needs of the user.

Sometimes this method is applied to match large and small numbers:

abs(X - Y) / (1 + min(abs(X), abs(Y)))

[/EDITED]

But for real problems, neither the relative nor the absolute error are accurate in general:

There is no unique mathematically meaningful definition of "near to eachother" for two values. E.g. 1000001 and 1000000 can be assumed to be near to eachother using a relative tolerance. But 0.000001 and 1.000001 have the same distance and deciding if they are near or far apart depends on the physical definition of the value. If it is the temperature of water, the different values mean a totally different behaviour. If the values are the height of a measurement point relative to the sea level, measured in millimeter from a satellite, both values are almost equal.

In real scientific projects the metric used for the measurement of the local discretization error is very important. The integrator adjustes the step size to keep the local error under, but near to the specified tolerance. This sounds trivial, but imagine the simulation of the re-entrance of a spaceship and the trajectory is this vector:

y = [position; velocity; temperature inside; temperature surface]

While the temperature are relevant for the application, they do not influence the spacial motion of the body. In consequence the local error should not be controlled neither by relative or absolute tolerances of the temperatures! In addition it should not matter, if the position is measured in AU or meter, or relative to the sea level or to the center of the earth.

Conclusion: For an accurate integration the definition of the local error must consider the physical meaning of the components of the trajectory. Using AbsTol, RelTol and NormControl (see help odeset) is not powerful enough for an accurate integration. A better method is either a user-defined function to calculate a norm, or local scaling, which let all relevant values live near to 1.0.

[Julián Francisco]

@Jan Simon: Thank you for your answer. I have forgotten to say the above Matlab definitions have been taken from help odeset. I waited an answer more focused on the questions I asked.

[Jan]

See EDITED section for some more details. I assume the definition gets clearer, when you try it with some example values.

[Julián Francisco]

@Jan Simon:

- I asked "why RelTol controls the number of correct digits in all solution components, except those smaller than thresholds AbsTol(i)". What happens if any of the solution components are smaller than its corresponding AbsTol value? Does AbsTol give the number of correct digits in this case?

- I am using the default integrators of Matlab (ode23, ode45 and ode113) for integrating my equations. What others integrators or methods to do it do you recommend me?

[Jan]

@Julian: If your trajectory has two components y=[y1; y2] and while y1 is near to 7000 (e.g. the position), y2 is near to 0 (e.g. the speed). Then controlling the local discretization error inside the integrator cannot use the same method for both components. Therefore for [y>AbsTol] the relative tolerance is used, and for [y<=AbsTol] the absolute tolerance. As you can see by interpreting the formulas or by a tiny numerical experiment, the absolute tolerance does *not necessarily* determine the correct number of digits. It does depend on the decision, what exactly "near" to zero means for your physical object: If you measure the height of a supertanker above the sea-level, 1mm and 0.00001mm are very near and the number of "valid" digits is determined by the absolute tolerance. But if you measure the thickness of a goldfilm, 1mm and 0.00001mm are totally different.

This means, that even then question is not logical. If the trajectory is "near" to zero, the term "correct digits" is not well-defined anymore.

I'm using a hand-coded ODE45 integrator, which let the user define a norm-function for the local discretization error and event functions, which can influence the parameters. In addition it can use the sensitivity matrix (dependency of the trajectory to small variations of the initial value) to control the stepsize. Unfortunately I do not find the time to polish the function for posting it in the FEX. Therefore aksing the technical support would be a good idea. If they are fast, an improved integrator can be included in Matlab 2012b.

[Julián Francisco]

@Jan Simon: Thank you very much for your comments.

Answer 2

[Mike Hosea]

The numerical method works with an error estimate, i.e. it computes an approximation for abs(x-x0), where x is the approximate solution and x0 is the exact solution. Let's call the error estimate E. The code will try to meet E < max(abstol,reltol*abs(x)). When abs(x) is larger than abstol/reltol, then only the relative error tolerance is important, so you are asking that E/abs(x) < reltol. E/abs(x) is an approximation of the relative error abs(x-x0)/abs(x0). So you essentially control the number of significant digits, as you have noted. However, when abs(x) is smaller than abstol/reltol, only the absolute error tolerance is used, so the error test in that case is E < abstol, which is approximately abs(x-x0) < abstol. Since relative error is undefined when the true solution is zero, you have to switch over to something besides relative error control when the solution component gets small. E.g. if x0 = 0, the approximation x = realmin has no correct significant digits, nary a one.

[Mike Hosea]

BTW, what is "small" depends on the scaling of the component. Quite generally, near zero results arise from cancellation, and the absolute tolerance that makes sense depends on the sizes of the numbers before cancellation. If a solution component is scaled small, e.g. around 1e-20, and stays in that range, then presumably you can set abstol(i) for that component to 1e-26 to get about 6 digits of accuracy. However, this will probably be impossible if the small value of this component arises because much larger numbers cancel.

[Jan]

@Mike: Exactly. See my example of the supertanker: If the height is measured from the center of earth, the relative error determines the correct number of digits. If it is measured from the sea-level (cancellation!), the absolute error is required. +1

[Julián Francisco]

@Mike Hosea: Thank you so much for your answer and comments.