Working with Complex Numbers in Julia

7.3. Working with Complex Numbers in Julia#

Complex numbers are a built-in numeric type in Julia. The global constant im is bound to the imaginary unit, \(\sqrt{-1}\).

7.3.1. Creating Complex Numbers#

You can create a complex number by combining a real part and an imaginary part. Note that the coefficient must directly multiply im.

# The imaginary unit `im` is a suffix to a numeric literal.
a = 5 - 3im

5 - 3im

Alternatively, you can use the complex() constructor function.

# This syntax is useful when the real and imaginary parts are variables.
b = complex(5, 3)

5 + 3im

When the imaginary part is a fraction or an expression, you need to use parentheses. Julia’s order of operations would otherwise interpret 3/4im as 3/(4*im).

# Parentheses are required to group the coefficient of `im`.
c = 5 + (3/4)im

5.0 + 0.75im

7.3.2. Standard Arithmetic Operations#

All standard arithmetic operations are defined for complex numbers and work as you’d expect from mathematics.

# (5 - 3i) * (5 + 3i) = 25 + 15i - 15i - 9i^2 = 25 + 9 = 34
a * b

34 + 0im

# Division is also supported.
a / b

0.4705882352941177 - 0.8823529411764706im

# (1 + i)^2 = 1 + 2i + i^2 = 1 + 2i - 1 = 2i
(1 + im)^2

0 + 2im

7.3.3. Manipulating Complex Numbers#

Julia provides a standard library of functions to work with the different parts of a complex number.

z = 1 + 2im

1 + 2im

# Get the real part of z
real(z)

1

# Get the imaginary part of z
imag(z)

2

# Get the complex conjugate of z
conj(z)

1 - 2im

# Get the absolute value (magnitude or modulus) of z
abs(z)

2.23606797749979

# Get the squared absolute value, which avoids a square root and is often more efficient.
abs2(z)

5

# Get the phase angle in radians.
angle(z)

1.1071487177940904

7.3.4. Elementary Functions with Complex Arguments#

Most elementary functions are also defined for complex arguments, enabling a wide range of mathematical computations.

sqrt(1im)

0.7071067811865476 + 0.7071067811865475im

cos(1 + 2im)

2.0327230070196656 - 3.0518977991517997im

exp(1 + 2im)

-1.1312043837568135 + 2.4717266720048188im

An important principle in Julia is type stability: functions typically return a value of a predictable type based on the input types. For this reason, sqrt will only return a complex number if its input is complex. Giving it a negative real number results in a DomainError.

# This is fine: complex input gives a complex output.
sqrt(-1.0 + 0im)

0.0 + 1.0im

# This will error: the input is a real number, so the output must be real.
# The error message helpfully suggests how to fix this.
sqrt(-1.0)

DomainError with -1.0:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).

Stacktrace:
 [1] throw_complex_domainerror(f::Symbol, x::Float64)
   @ Base.Math ./math.jl:33
 [2] sqrt(x::Float64)
   @ Base.Math ./math.jl:608
 [3] top-level scope
   @ In[18]:3

7.3.5. Application: A Robust Quadratic Formula#

We can now write a single, robust function to find the roots of a quadratic equation, \(ax^2+bx+c=0\). By ensuring the discriminant (\(b^2 - 4ac\)) is converted to a complex number before taking the square root, our function will seamlessly handle both real and complex roots.

function roots_of_quadratic(a, b, c)
    # Calculate the discriminant.
    discriminant = b^2 - 4*a*c
    
    # Convert the discriminant to a complex number before taking the square root.
    # This is the key step that handles all cases!
    d = sqrt(complex(discriminant))
    
    r1 = (-b - d) / (2a)
    r2 = (-b + d) / (2a)
    
    return r1, r2
end

roots_of_quadratic (generic function with 1 method)

Let’s test it on an equation with complex roots, \(x^2 + 1 = 0\).

roots_of_quadratic(1, 0, 1)

(0.0 - 1.0im, 0.0 + 1.0im)

And another one, \(x^2 + 4x + 13 = 0\).

roots_of_quadratic(1, 4, 13)

(-2.0 - 3.0im, -2.0 + 3.0im)

Because we always convert the discriminant to a complex number, our function will return complex numbers even if the roots are real. Notice the imaginary part is 0.0im.

# Testing with -x^2 + 5x + 6 = 0, which has real roots.
roots_of_quadratic(-1, 5, 6)

(6.0 + 0.0im, -1.0 - 0.0im)

This is often the desired behavior for consistency. If you needed to, you could modify the function to check if the imaginary part is zero and return real numbers in that case.

Contents