Complex

// No real number squares to -1, so `x * x + 1 = 0` has no real solution. Complex numbers add
// one: `i`, with `i * i = -1`. A complex number is a pair, `real + imaginary·i`, and with it
// every quadratic has its two roots, and a rotation in the plane becomes a multiplication.
//
// The Math package's `Complex` holds the pair as two `float64` fields, `real` and `imaginary`,
// and gives it the four operators, so complex arithmetic is written as ordinary arithmetic.
// (`Complex32` is the same over `float32`.)
//
// The pair can also be read as a point in the plane, or as an arrow from the origin to it. The
// arrow's length is the magnitude, `Abs()`, and its direction is the angle, `Arg()`, measured in
// radians from the positive real axis. `FromPolar` builds a number from those two instead.
import Io::PrintLine;
import Math::{ Complex, DegToRad, Pi, RadToDeg };

// A Complex does not implement Display, so `PrintLine("{}", a)` is refused:
//
//     error: no matching overload for 'PrintLine' with argument types (char8[..], Complex)
//
// This writes it as `a + bi` or `a - bi` instead.
func Show(label: char8[..], z: Complex) {
    if z.imaginary < 0.0 {
        PrintLine("{:<14} {:.3} - {:.3}i", label, z.real, -z.imaginary);
    } else {
        PrintLine("{:<14} {:.3} + {:.3}i", label, z.real, z.imaginary);
    }
}

func Main() -> int {
    // Arithmetic. Multiplying uses i·i = -1: (3 + 4i)(1 - 2i) = 3 - 6i + 4i + 8 = 11 - 2i.
    let a = Complex::New(3.0, 4.0);
    let b = Complex::New(1.0, -2.0);
    Show("a + b", a + b);
    Show("a * b", a * b);
    Show("a / b", a / b);
    Show("conjugate a", a.Conjugate());

    // The operators take a Complex on both sides; a plain float is refused:
    //
    //     error: operator '*' cannot combine left operand 'Complex' with right operand 'float64'
    //
    // `Scale` multiplies by a real number instead.
    Show("a scaled by 2", a.Scale(2.0));

    // Polar form: 3 + 4i is 5 long and points about 53 degrees above the real axis.
    PrintLine("{:<14} {:.3}", "magnitude", a.Abs());
    PrintLine("{:<14} {:.3} degrees", "angle", RadToDeg(a.Arg()));

    // Multiplying by i turns a number a quarter turn anticlockwise, keeping its length.
    let i = Complex::New(0.0, 1.0);
    let turned = a * i;
    Show("a * i", turned);
    PrintLine("{:<14} {:.3}, {:.3} degrees", "turned polar", turned.Abs(), RadToDeg(turned.Arg()));

    // `FromPolar` goes the other way: length 2 at 60 degrees.
    Show("polar 2, 60", Complex::FromPolar(2.0, DegToRad(60.0)));

    // x^2 + 2x + 5 = 0. The discriminant, 4 - 20, is negative, so the real square root would be
    // NaN; the complex one is 4i, and the two roots are -1 + 2i and -1 - 2i.
    let p = 1.0;
    let q = 2.0;
    let r = 5.0;
    let root = Complex::New(q * q - 4.0 * p * r, 0.0).Sqrt();
    let first = (Complex::New(-q, 0.0) + root).Scale(1.0 / (2.0 * p));
    let second = (Complex::New(-q, 0.0) - root).Scale(1.0 / (2.0 * p));
    Show("first root", first);
    Show("second root", second);

    // Putting a root back in gives zero, as far as floats allow.
    let check = first * first + first.Scale(q) + Complex::New(r, 0.0);
    Show("check", check);

    // Euler's identity, e^(i·pi) = -1. The imaginary part comes out as about 1.2e-16, not 0,
    // because `Pi` is only the nearest float to pi; at three places it shows as 0.000.
    let euler = Complex::New(0.0, Pi).Exp();
    Show("e^(i pi)", euler);
    PrintLine("{:<14} {}", "imaginary part", euler.imaginary);
    return 0;
}