$$ \newcommand{\tp}{\thinspace .} $$

This chapter is taken from the book A Primer on Scientific Programming with Python by H. P. Langtangen, 5th edition, Springer, 2016.

Complex numbers

Suppose \( x^2=2 \). Then most of us are able to find out that \( x=\sqrt{2} \) is a solution to the equation. The more mathematically interested reader will also remark that \( x=-\sqrt{2} \) is another solution. But faced with the equation \( x^2=-2 \), very few are able to find a proper solution without any previous knowledge of complex numbers. Such numbers have many applications in science, and it is therefore important to be able to use such numbers in our programs.

On the following pages we extend the previous material on computing with real numbers to complex numbers. The text is optional, and readers without knowledge of complex numbers can safely drop this section and jump to the section Summary.

A complex number is a pair of real numbers \( a \) and \( b \), most often written as \( a+bi \), or \( a+ib \), where \( i \) is called the imaginary unit and acts as a label for the second term. Mathematically, \( i=\sqrt{-1} \). An important feature of complex numbers is definitely the ability to compute square roots of negative numbers. For example, \( \sqrt{-2}=\sqrt{2}i \) (i.e., \( \sqrt{2}\sqrt{-1} \)). The solutions of \( x^2=-2 \) are thus \( x_1=+\sqrt{2}i \) and \( x_2=-\sqrt{2}i \).

There are rules for addition, subtraction, multiplication, and division between two complex numbers. There are also rules for raising a complex number to a real power, as well as rules for computing \( \sin z \), \( \cos z \), \( \tan z \), \( e^z \), \( \ln z \), \( \sinh z \), \( \cosh z \), \( \tanh z \), etc. for a complex number \( z = a+ib \). We assume in the following that you are familiar with the mathematics of complex numbers, at least to the degree encountered in the program examples. $$ \begin{equation*} \hbox{let } u=a+bi\hbox{ and }v=c+di\end{equation*} $$

The following rules reflect complex arithmetics: $$ \begin{align*} u &= v\quad\Rightarrow\quad a=c, \ b=d\\ -u &= -a -bi\\ u^* &\equiv a - bi\quad\hbox{(complex conjugate)}\\ u+ v &= (a+c) + (b+d)i\\ u - v &= (a-c) + (b-d)i\\ uv &= (ac - bd) + (bc+ad)i\\ u/v &= \frac{ac +bd}{c^2 + d^2} + \frac{bc-ad}{c^2+d^2}i\\ |u| &= \sqrt{a^2 + b^2}\\ e^{iq} &= \cos q + i\sin q \end{align*} $$

Complex arithmetics in Python

Python supports computation with complex numbers. The imaginary unit is written as j in Python, instead of \( i \) as in mathematics. A complex number \( 2-3i \) is therefore expressed as (2-3j) in Python. We remark that the number \( i \) is written as 1j, not just j. Below is a sample session involving definition of complex numbers and some simple arithmetics:

>>> u = 2.5 + 3j       # create a complex number
>>> v = 2              # this is an int
>>> w = u + v          # complex + int
>>> w
(4.5+3j)

>>> a = -2
>>> b = 0.5
>>> s = a + b*1j       # create a complex number from two floats
>>> s = complex(a, b)  # alternative creation
>>> s
(-2+0.5j)
>>> s*w                # complex*complex
(-10.5-3.75j)
>>> s/w                # complex/complex
(-0.25641025641025639+0.28205128205128205j)
A complex object s has functionality for extracting the real and imaginary parts as well as computing the complex conjugate:

>>> s.real
-2.0
>>> s.imag
0.5
>>> s.conjugate()
(-2-0.5j)

Complex functions in Python

Taking the sine of a complex number does not work:

>>> from math import sin
>>> r = sin(w)
Traceback (most recent call last):
  File "<input>", line 1, in ?
TypeError: can't convert complex to float; use abs(z)
The reason is that the sin function from the math module only works with real (float) arguments, not complex. A similar module, cmath, defines functions that take a complex number as argument and return a complex number as result. As an example of using the cmath module, we can demonstrate that the relation \( \sin(ai)=i\sinh a \) holds:

>>> from cmath import sin, sinh
>>> r1 = sin(8j)
>>> r1
1490.4788257895502j
>>> r2 = 1j*sinh(8)
>>> r2
1490.4788257895502j
Another relation, \( e^{iq} = \cos q + i\sin q \), is exemplified next:

>>> q = 8     # some arbitrary number
>>> exp(1j*q)
(-0.14550003380861354+0.98935824662338179j)
>>> cos(q) + 1j*sin(q)
(-0.14550003380861354+0.98935824662338179j)

Unified treatment of complex and real functions

The cmath functions always return complex numbers. It would be nice to have functions that return a float object if the result is a real number and a complex object if the result is a complex number. The Numerical Python package has such versions of the basic mathematical functions known from math and cmath. By taking a

from numpy.lib.scimath import *
one obtains access to these flexible versions of mathematical functions. The functions also get imported by any of the statements

from scipy import *
from scitools.std import *

A session will illustrate what we obtain. Let us first use the sqrt function in the math module:

>>> from math import sqrt
>>> sqrt(4)     # float
2.0
>>> sqrt(-1)    # illegal
Traceback (most recent call last):
  File "<input>", line 1, in ?
ValueError: math domain error
If we now import sqrt from cmath,

>>> from cmath import sqrt
the previous sqrt function is overwritten by the new one. More precisely, the name sqrt was previously bound to a function sqrt from the math module, but is now bound to another function sqrt from the cmath module. In this case, any square root results in a complex object:

>>> sqrt(4)     # complex
(2+0j)
>>> sqrt(-1)    # complex
1j
If we now take

>>> from numpy.lib.scimath import *
we import (among other things) a new sqrt function. This function is slower than the versions from math and cmath, but it has more flexibility since the returned object is float if that is mathematically possible, otherwise a complex is returned:

>>> sqrt(4)     # float
2.0
>>> sqrt(-1)    # complex
1j

As a further illustration of the need for flexible treatment of both complex and real numbers, we may code the formulas for the roots of a quadratic function \( f(x)=ax^2 + bx + c \):

>>> a = 1; b = 2; c = 100   # polynomial coefficients
>>> from numpy.lib.scimath import sqrt
>>> r1 = (-b + sqrt(b**2 - 4*a*c))/(2*a)
>>> r2 = (-b - sqrt(b**2 - 4*a*c))/(2*a)
>>> r1
(-1+9.94987437107j)
>>> r2
(-1-9.94987437107j)
Using the up arrow, we may go back to the definitions of the coefficients and change them so the roots become real numbers:

>>> a = 1; b = 4; c = 1   # polynomial coefficients
Going back to the computations of r1 and r2 and performing them again, we get

>>> r1
-0.267949192431
>>> r2
-3.73205080757
That is, the two results are float objects. Had we applied sqrt from cmath, r1 and r2 would always be complex objects, while sqrt from the math module would not handle the first (complex) case.