Proposed algorithm

We an algorithm for the reciprocal scaling of a complex vector X by a complex number A. We use the notation AR and AI to denote the real and imaginary parts of A, respectively. The algorithm is as follows:

1. If AI is zero, then call {CS,ZD}RSCL (currently in LAPACK) using AR.
2. If AR is zero, then if AI is in the safe range, call SCAL with the number CMPLX( ZERO, -ONE / AI ). Otherwise, do the proper scaling by a power of two similarly to CSRSCL.
3. If both the real and imaginary parts of A are nonzero, then we compute UR and UI such that CMPLX( ONE / UR, - ONE / UI ) is the reciprocal of A, i.e.,

   UR = AR + AI * ( AI / AR )
   UI = AI + AR * ( AR / AI )

   Note that UR and UI are always different from zero. NaNs only appear if either:
   1. AR or AI is a NaN.
   2. AR and AI are both infinite, in which case it makes sense to propagate a NaN.
      Now, we have three cases:
   3. If UR and UI are both in the safe range, then call SCAL with the number CMPLX( ONE / UR, - ONE / UI ).
   4. If the absolute values of either UR or UI are smaller than SFMIN, it means that both UR and UI are small numbers. In fact, we can prove that they should be both smaller than SFMIN(1+1/EPSILON) Therefore, it makes sense to scale x by CMPLX( SFMIN / UR, - SFMIN / UI ), and then scale the result by ONE / SFMIN. Use SCAL to do both scalings.
   5. If the absolute values of either UR or UI are greater than ONE / SFMIN, then we must check a few things:
      1. If either AR or AI is infinite, then UR and UI are either both infinite or both NaNs. Therefore, we scale by CMPLX( ONE / UR, - ONE / UI ) using SCAL.
      2. If either UR or UI is infinite and AR and AI are finite, then the algorithm generated infinite numbers that can be avoided. In this case, recompute scaled versions of UR and UI using SFMIN. Then, scale x by SFMIN, and then scale the result by CMPLX( ONE / UR, - ONE / UI ). Use SCAL to do both scalings.
      3. In other cases, scale x by SFMIN, and then scale the result by CMPLX( ONE / (SFMIN*UR), - ONE / (SFMIN*UI) ). Use SCAL to do both scalings.

Some comments about the algorithm:

• It does not use complex divisions. Instead, it uses $n$ complex multiplications, $n$ real multiplications, and a fixed number of real multiplications and divisions.
• It propagate NaNs in AR and AI to the result. Moreover, it propagates NaNs if AR and AI are both infinite. It does not generate NaNs in other cases.
• It does not generate infinity numbers unless A is zero.
• It avoids generating zeros in both the real and imaginary parts of the reciprocal of A. This is done by treating cases where UR or UI are infinite numbers.
• It does not lose accuracy due to subnormal numbers when generating the reciprocal of A.