Texas Instrument CLA float inverse trick - What is its purpose

Viewed 104

I have this piece of code written by somebody else that runs on a TI TMS320 Command Law Accelerator. So it's optimized in size and speed.

In order to get 1/x, the code always does something like this.

float32 y = __meinvf32(x); 
y = y * (2.0f - y*x); 
y = y * (2.0f - y*x);

I found this thread that propose something similar, but in my case, there is no clamping at the end.

Can seomebody help me understand what is the intent behind this?

2 Answers

Isaac Newton figured this out.

__meinvf32(x) gives an approximation of 1/x, say 1/x • (1+e), where e is some small relative error.

Let y = 1/x • (1+e). Then, when we calculate y • (2 − yx), we have:

  • y • (2 − yx) =
  • (1/x • (1+e)) • (2 − (1/x • (1+e))•x) =
  • 1/x • (1+e) • (2 − (1+e)) =
  • 1/x • (2 + 2e − (1+e) − e(1+e)) =
  • 1/x • (2 + 2e − 1 − eee2) =
  • 1/x • (1 − e2).

Since e is small, e2 is even smaller. Thus, by calculating y • (2 − yx) we get an estimate of 1/x that is closer than before; the relative error is only −e2 instead of e. Repeating this improves the estimate again (up to the limits of floating-point precision).

With some knowledge of bounds on the initial e, we can calculate how many repetitions are needed to get the estimate as close as desired to the correct result.

y = e + 1/x where e is some small error.

So (2.0 - y*x) is close to 1.0 and has the effect of reducing e on each pass.

Related