# How to compute the square root of a complex number

Suppose you’re given a complex number

z = x + iy

and you want to find a complex number

w = u + iv

such that w² = z. If all goes well, you can compute w as follows:

ℓ = √(x² + y²)
u = √((ℓ + x)/2)
v = sign(y) √((ℓ – x)/2).

For example, if z = 5 + 12i, then ℓ = 13, u = 3, and v = 2. A quick calculation confirms

(3 + 2i)² = 5 + 12i.

(That example worked out very nicely. More on why in the next post.)

Note that ℓ ≥ |x|, and so the square roots above have non-negative arguments.

## Numerical issues

I said “if all goes well” above because in floating point arithmetic, various things could go wrong. In , the authors point out that x² + y² could overflow. They give a more complicated algorithm that addresses the problem of overflow, and they reference a paper with an even more complicated algorithm to cover all the corner cases — NaNs, infinities, etc.

If you have access to a `hypot` function, it will let you compute √(x² + y²) without overflow. If you don’t have access to a `hypot` function, you can roll your own quickly using the algorithm here.

This may seem unnecessary. Why would you need to worry about overflow? And why would you have access to a real square root function but not a complex square root function?

If you’re using standard 64-bit floating point numbers, the hypotenuse calculation is not going to overflow as long as your arguments are less than 10154, and you can use the function `csqrt` from `complex.h`.

But if you’re using a 16-bit float or even an 8-bit float, overflow is an ever-present concern, and you might not have the math library support you’re accustomed to.

## Related posts

 Jean-Michel Muller et al. Handbook of Floating Point Arithmetic.