f(z) Complex Numbers and Functions Tech Note Complex Numbers

Humans invented several number systems through the ages, including

1. Set of integers constructed from the counting numbers, 1, 2, 3, …
2. Rational numbers, m/n, as ratio of integers.
3. Real numbers, x, including rational numbers (e.g., ½) and irrational numbers (e.g., ).
4. Complex numbers, a + ib, where a is the “real part” and ib is the “imaginary part.”

Every “new” set of invented numbers extended the previous system.

These web pages are a brief overview of complex numbers and complex math, and describe a ComplexMathLibrary unit that has a number of functions useful for numerical experiments.  Unlike many other complex math libraries, many algorithms in ComplexMathLibrary are implemented in both rectangular and polar forms, since some algorithms are easier to implement in one form than the other.

This page gives an introduction to imaginary and complex numbers both in rectangular and polar form, and an overview of the software implementation.  Other pages explain complex arithmetic, complex logarithms and powers, complex trig functions, complex hyperbolic functions, complex Bessel functions, and the complex gamma function.

“Imaginary” Numbers

The solution to the equation, x2 = 1, is

But what is the solution to the equation, x2 = -1?  By definition, we call the “imaginary” number “i” the solution to this equation:

(Electrical engineers often use a “j” instead of an “i” to avoid confusion with the abbreviation “i” being used as “electric current” in equations such as Ohm’s Law.)

Imaginary numbers are “different” from real numbers and cannot be combined with real numbers.  So, 1 + 2 = 3, and i + 2i = 3i.  But 1 + 2i is something else, namely a “complex” number.

Complex Numbers

An ordered pair, (x, y), of a real number, x, and an imaginary number, iy, is called a “complex” number.  A complex number, z, can be written as z = (x, y), or z = x + iy, (or in other forms, which will be discussed below).

The complex plane is a two-dimensional space, where each point can be identified in (x, y) rectangular coordinates, or (r, q) polar coordinates.

Rectangular Form

The absolute value (or sometimes called modulus or “norm”) of the complex value z is defined as follows:

Results from CAbs(z) function

 Absolute value of complex number:  CAbs = ABS(z)      z                ABS(z) ------------  -------------------  0.0 +  0.0i    0.000000000000000   0.5 +  0.5i    0.707106781186548  -0.5 +  0.5i    0.707106781186548  -0.5 -  0.5i    0.707106781186548   0.5 -  0.5i    0.707106781186548   1.0 +  0.0i    1.000000000000000   1.0 +  1.0i    1.414213562373095   0.0 +  1.0i    1.000000000000000  -1.0 +  1.0i    1.414213562373095  -1.0 +  0.0i    1.000000000000000  -1.0 -  1.0i    1.414213562373095   0.0 -  1.0i    1.000000000000000   1.0 -  1.0i    1.414213562373095   5.0 +  0.0i    5.000000000000000   5.0 +  3.0i    5.830951894845301   0.0 +  3.0i    3.000000000000000  -5.0 +  3.0i    5.830951894845301  -5.0 +  0.0i    5.000000000000000  -5.0 -  3.0i    5.830951894845301   0.0 -  3.0i    3.000000000000000  -5.0 -  3.0i    5.830951894845301

Polar Form

Instead of (x,y) rectangular coordinates, all points can be represented in polar form with (r, q) coordinates as shown in the diagram below:

Normally the angle is limited to a symmetric range, -p < q £ p, instead of 0 < q £ 2p.

Sometimes the form, z = r[cos q + i·sin q], is abbreviated as z = r·cis q.

Software Implementation

There are several approaches that can be used in implementing a ComplexMathLibrary unit for Delphi/Kylix.  I view a complex data type much like I would an integer or real numeric data type:  a low-level primitive type.  Some may choose to use an “object oriented” approach and implement a TComplex data type as an “object.”   I prefer to use such an object-oriented approach only for a “macro” object (perhaps a matrix object), and use a more primitive “Record” to define a TComplex type, since from a mathematical standpoint a complex value “behaves” much like a primitive floating-point value.

Some complex math algorithms can be expressed more easily using polar coordinates instead of rectangular coordinates.  The TComplex data type shown below allows definition of a complex value in either rectangular or polar coordinates.  Conversions between the types are forced when necessary (like to add two complex numbers in polar form), but otherwise the programmer can control what format is appropriate for a given algorithm.  The “Case” statement in the “Record” definition allows either (x, y) rectangular coordinates, or (r, q) polar coordinates to be stored in the same data structure.

TComplex data type

 TYPE   TReal        = DOUBLE;   TComplexForm = (cfPolar, cfRectangular);     TComplex     =     RECORD       CASE form:  TComplexForm OF         cfRectangular:  (x,y    :  TReal);  // z = x + i*y         cfPolar      :  (r,theta:  TReal);  // z = r*CIS(theta)     END;                // where CIS(theta) = COS(theta) + i*SIN(theta)                         //       theta = -PI..PI (in canonical form)

The assignment of a TComplex variable is usually via a CSet function:

CSet function

 FUNCTION CSet (CONST a,b:  TReal; CONST f:  TComplexForm):  TComplex; BEGIN   RESULT.form := f;   CASE f OF     cfRectangular:       BEGIN         RESULT.x := a;         RESULT.y := b       END;     cfPolar:       BEGIN         RESULT.r := a;         RESULT.theta := b       END   END END {CSet};

The function prototype for CSet defines a default TComplexForm of cfRectangular.

The array, a, of complex values used in many of the ComplexMathVerification project routines was defined as follows:

Definition of complex array

 VAR a:  ARRAY[1..21] OF TComplex; … a[ 1] := CSet(  0.0,  0.0); a[ 3] := CSet( -0.5,  0.5); a[ 4] := CSet( -0.5, -0.5); a[ 5] := CSet(  0.5, -0.5); a[ 6] := CSet(  1.0,  0.0); a[ 7] := CSet(  1.0,  1.0); a[ 8] := CSet(  0.0,  1.0); a[ 9] := CSet( -1.0,  1.0); a[12] := CSet(  0.0, -1.0); a[13] := CSet(  1.0, -1.0); a[16] := CSet(  0.0,  3.0); a[17] := CSet( -5.0,  3.0); a[19] := CSet( -5.0, -3.0); a[20] := CSet(  0.0, -3.0); a[21] := CSet( -5.0, -3.0);   // "Polar" definitions a[ 2] := CSet(SQRT(0.5),       PI/4, cfPolar);   //  0.5 + 0.5i a[10] := CSet(     1.0,          PI, cfPolar);   // -1   + 0i a[11] := CSet( SQRT(2),      5*PI/4, cfPolar);   // -1   -  i a[14] := CSet(     5.0,         0.0, cfPolar);   //  5   + 0i a[15] := CSet(SQRT(34), ArcTan(3/5), cfPolar);   //  5   + 3i a[18] := CSet(     5.0,          PI, cfPolar);   // -5   + 0i

Such complex values can be displayed in either rectangular or polar form using the CToRectStr and CToPolarStr functions (explained below).

Code to display array of complex values using CToRectStr and CToPolarStr

 Memo.Lines.Add('Complex number definition and conversion:  ' +                      'CSet, CtoRectStr, CtoPolarStr'); Memo.Lines.Add(''); Memo.Lines.Add('             rectangular          ' +                '              polar'); Memo.Lines.Add('     ----------------------------   ' +                '--------------------------------'); FOR k := Low(a) TO High(a) DO BEGIN   Memo.Lines.Add(   Format('%2d  %-s %-s',                             [k, CToRectStr(a[k],  13,9),                                 CToPolarStr(a[k], 13,9)])) END;

The code above results in the follow TMemo data:

Array of TComplex values

 Complex number definition and conversion:  CSet, CtoRectStr, CtoPolarStr                rectangular                        polar      ----------------------------   --------------------------------  1    0.000000000 +  0.000000000i   0.000000000 * CIS(  0.000000000)  2    0.500000000 +  0.500000000i   0.707106781 * CIS(  0.785398163)  3   -0.500000000 +  0.500000000i   0.707106781 * CIS(  2.356194490)  4   -0.500000000 -  0.500000000i   0.707106781 * CIS( -2.356194490)  5    0.500000000 -  0.500000000i   0.707106781 * CIS( -0.785398163)  6    1.000000000 +  0.000000000i   1.000000000 * CIS(  0.000000000)  7    1.000000000 +  1.000000000i   1.414213562 * CIS(  0.785398163)  8    0.000000000 +  1.000000000i   1.000000000 * CIS(  1.570796327)  9   -1.000000000 +  1.000000000i   1.414213562 * CIS(  2.356194490) 10   -1.000000000 +  0.000000000i   1.000000000 * CIS(  3.141592654) 11   -1.000000000 -  1.000000000i   1.414213562 * CIS(  3.926990817) 12    0.000000000 -  1.000000000i   1.000000000 * CIS( -1.570796327) 13    1.000000000 -  1.000000000i   1.414213562 * CIS( -0.785398163) 14    5.000000000 +  0.000000000i   5.000000000 * CIS(  0.000000000) 15    5.000000000 +  3.000000000i   5.830951895 * CIS(  0.540419500) 16    0.000000000 +  3.000000000i   3.000000000 * CIS(  1.570796327) 17   -5.000000000 +  3.000000000i   5.830951895 * CIS(  2.601173153) 18   -5.000000000 +  0.000000000i   5.000000000 * CIS(  3.141592654) 19   -5.000000000 -  3.000000000i   5.830951895 * CIS( -2.601173153) 20    0.000000000 -  3.000000000i   3.000000000 * CIS( -1.570796327) 21   -5.000000000 -  3.000000000i   5.830951895 * CIS( -2.601173153)

A dynamic array of TComplex values works just like any other Delphi (D4 or later) dynamic array:

CSet with 2D dynamic array of TComplex values

 VAR  z:  ARRAY OF ARRAY OF TComplex; . . . BEGIN   SetLength(z, ImageComposite.Height, ImageComposite.Width);   FOR j := Low(z) TO High(z) DO   BEGIN     y := 3 - 6*(j-Low(z))/(High(z) - Low(z));     FOR i := Low(z[j]) TO High(z[j]) DO     BEGIN       x := -4 + 8*(i-Low(z[j]))/(High(z[j]) - Low(z[j]));       z[j,i] := CSet(x,y, cfRectangular);     END   END; . . .

Converting rectangular coordinates to polar coordinates

Conversion from (x, y) rectangular form to polar coordinates, (r, q), is fairly simple:

To provide a range of q  to be in the interval [-p, p), the arctan2(y, x) function is used from the Delphi math library.

Converting polar coordinates to rectangular coordinates

The opposite conversion, from (r, q) polar coordinates back to the (x, y) rectangular coordinates, is also simple:

The CConvert function can be used to convert from one format to the other:

CConvert function

 // complex number definition/conversion/output FUNCTION CConvert (CONST z:  TComplex; f:  TComplexForm):  TComplex;   VAR     a:  TComplex; BEGIN   IF   IsInfinity(z.x) OR IsInfinity(z.y)   THEN BEGIN     CASE f OF       cfPolar:        RESULT := CSet(Infinity, 0.0, cfPolar);       cfRectangular:  RESULT := CSet(Infinity, Infinity)     END   END   ELSE BEGIN     IF   IsNan(z.x) OR IsNan(z.y)     THEN RESULT := CSet(NaN, NaN)     ELSE BEGIN         IF   z.form = f       THEN RESULT := z       ELSE BEGIN         CASE z.form OF           cfPolar:            // cfPolar-to-cfRectangular conversion             BEGIN               a.form := cfRectangular;               a.x := z.r * COS(z.theta);               a.y := z.r * SIN(z.theta)             END;             cfRectangular:      // cfRectangular-to-cfPolar conversion             BEGIN               a.form := cfPolar;               a.r := Cabs(z);               a.theta := ARCTAN2(z.y, z.x);             END;         END;         RESULT := a       END     END   END END {CConvert};

Using CConvert, any TComplex value can be printed either rectangular or polar format.  The default format for either form is specified in the FormatString in the function header:

CToPolarStr prototype with default string '%*.*f * CIS(%*.*f)'

 FUNCTION  CToPolarStr (CONST z:  TComplex; CONST w,d:  BYTE;                        CONST FormatString:  STRING =                         '%*.*f * CIS(%*.*f)'):  STRING;

CToPolar function

 FUNCTION  CToPolarStr (CONST z:  TComplex;                        CONST w,d:  BYTE;                        CONST FormatString:  STRING):  STRING;   VAR     u :  TComplex; BEGIN   u := CConvert (z, cfPolar);   RESULT := Format(FormatString, [w,d,u.r, w,d,u.theta]); END {CToPolarStr};

CToRectStr prototype with default string '%*.*f +%*.*fi'

 FUNCTION  CToRectStr  (CONST z:  TComplex; CONST w,d:  BYTE;                        CONST FormatString:  STRING =                         '%*.*f +%*.*fi'):  STRING;

CToRectStr function

 FUNCTION  CToRectStr (CONST z:  TComplex;                       CONST w,d:  BYTE;                       CONST FormatString:  STRING):  STRING;   VAR     NegativeFormat:  STRING;     u             :  TComplex; BEGIN   IF   IsNan(z.x) OR IsNan(z.y)   THEN RESULT := Format(FormatString, [w,d,NaN, w,d,NaN])   ELSE BEGIN     u := CConvert (z, cfRectangular);     IF   u.y >= 0     THEN RESULT := Format(FormatString, [w,d,u.x, w,d,u.y])     ELSE BEGIN       NegativeFormat := FormatString;       IF   POS('+', NegativeFormat) > 0       THEN NegativeFormat[ POS('+',NegativeFormat ) ] := '-';       RESULT := Format(NegativeFormat, [w,d,u.x, w,d,ABS(u.y)])     END   END END {CToRectStr};

Some extra work is made in CToRectStr so that two consecutive negative signs do not appear.

Euler’s formula and exponential form

Euler’s formula states:

One way to demonstrate the correctness of this formula is by using a series expansion of eiq, which breaks into real and imaginary series that are the cos and sin functions.

Only values on a unit circle can be described directly with Euler’s formula. Multiplication by some distance value, r, allows any complex value, z, to be described in this exponential form.  All of the following forms are equivalent:

Polar coordinates, (r, q), can be expressed in either a trigonometric or exponential form.