# Complex Function Library

1. Nov 6, 2012

### Enthalpy

Hello everybody!

Elsewhere someone searched in September a topic for her master thesis in software, and I suggested one in November , so I suppose the topic is available to other people...

It might be possible - perhaps - to improve the math function library for complex numbers.

The libraries I've seen add their layer to functions libraries on real numbers. This can be very inefficient: see how complicated it is to prevent all undue overflow errors in the complex division function given in "numerical receipes in C", a full page of tests and divisions. It is also indirect.

For the functions whose real version is not hard-wired but programmed or micro-programmed, maybe the algorithm can be adapted directly to complex numbers.

-----

Example: the real division is a multiplication by the reciprocal which is computed by an iterative algorithm, of numerical analysis type, often Newton-Raphson
http://en.wikipedia.org/wiki/Division_algorithm#Newton.E2.80.93Raphson_division
and this algorithm may well fit complex numbers with little effort.

-----

Example: the real square root is an iterative algorithm of numerical analysis type, like these
http://en.wikipedia.org/wiki/Methods_of_computing_square_roots#A_two-variable_iterative_method
http://en.wikipedia.org/wiki/Method...Iterative_methods_for_reciprocal_square_roots
again, such a method could run on complex numbers. Here Wiki seems to suggest the standard method
http://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Negative_or_complex_square
which would be pretty inefficient as the moduli demand already a square root.

-----

Example: the logarithm.

• It is computed on real numbers in base 2 using the floating point processor's normalisation. The exponent gives the heaviest bits and is discarded, then the mantissa is repeatedly squared and the exponent gives more bits and is discarded.
• Could it work on complex numbers? Harvest the biggest exponent of real or imaginary parts, subtract it from both. Iterate: square the complex, again harvest the biggest of both exponents, subtract from both.
• On complex numbers, I vaguely suppose that the bits of the next iteration's harvested exponent can overlap with the ones already obtained, but if you add the (shifted) contributions instead of concatenating them it may be exact and with no time penalty.
• By this way you obtain the logarithm of the modulus, which is the real part of the logarithm.

To get the imaginary part of the complex logarithm, which is the argument of the input complex, one can take advantage of the same successive squarings used to compute the real part.

When you square a complex, its argument (its angle, in polar coordinates) doubles, but modulo one turn. So on successive squarings, the sign of the imaginary part gives one bit of the argument at each step. You get the argument in units of turns.

Usable for atan2(r,i) as well but probably less interesting if you don't need a log.

-----

There may be other candidate functions: trigo, hyperbolic, exponential... More complicated and less frequent, most special functions are complex by nature.

I suppose the Newton-Raphson methods would transpose directly, but for a library, you have to check if they can underflow, overflow, how to minimise rounding errors... Maybe a mere normalisation at the beginning protects the whole function, which would again save time against a complex layer on real functions. Beware I didn't review all existing libraries.

Hardware (hi there) able to normalize two floats according to the biggest of both exponents (here considered the real and imaginary parts) would speed up the functions. It should bring the biggest part between 0.5 and 1, scale the other part identically, and return the exponent of the scaling factor. No idea if it's already available.

Also interesting, but I suppose it exists: an FPU instruction or micro-instruction that tells the quadrant, or better the octant, of a complex number, and brings the complex to the first quadrant or octant.

Marc Schaefer, aka Enthalpy

Share this great discussion with others via Reddit, Google+, Twitter, or Facebook

Can you offer guidance or do you also need help?
Draft saved Draft deleted