The minor nit first: the original specification did not support complex numbers, so adding them here is even less relevant than my return codes.

More importantly, your implementation is numerically less correct than my version, which itself is less correct than it should be.

Consider when b is near sqrt(d). In that case, b-sqrt(d) can lose precision. That's why I used the copysign function, so that I'm always adding two numbers of equal sign and size.

Mine is incorrect if b2 is near the size of 4ac since there too b2-4ac loses precision. Ideally this intermediate result should be done in quad precision if the input is in double.

Question: How does C++ handle copysign() (vs. copysignf() for floats), and support templates which want to use an higher precision intermediate value?