Implementing Gaussian integers in Haskell

Question

Currently my code (at GitHub) is this:

module Gaussian where

import qualified Data.Complex as Complex
  -- Enables toComplex.
  -- Qualified to avoid clashing.

data Gaussian = (:+)
  {
    real :: Integer,
    imaginary :: Integer
  } deriving
    (Show, Read, Eq)
  -- The gaussian integer type.
  -- A gaussian integer is a complex number of the form a + bi where a and b are integers.

instance Num Gaussian
  where
    (a1 :+ b1) + (a2 :+ b2) = (a1 + a2) :+ (b1 + b2)
    (a1 :+ b1) - (a2 :+ b2) = (a1 - a2) :+ (b1 - b2)
    (a1 :+ b1) * (a2 :+ b2) = (a1 * a2 - b1 * b2) :+ (a1 * b2 + a2 * b1)
    abs (a :+ b) = ((a ^ two) + (b ^ two)) :+ 0
      where
        two = 2 :: Integer
          -- This prevents compiler warnings about defaulting types.
    signum = id
    negate (a :+ b) = negate a :+ negate b
    fromInteger x = x :+ 0
  -- The Num instance for Gaussian.
  -- We use the euclidean norm as abs. norm(a + bi) = a^2 + b^2

i :: Gaussian
i = 0 :+ 1
  -- The imaginary unit, equal to sqrt(-1).

magnitude :: Floating c => Gaussian -> c
magnitude (a :+ b) = sqrt $ (fromInteger a ^ two) + (fromInteger b ^ two)
  where
    two = 2 :: Integer
      -- This prevents compiler warnings about defaulting types.

norm :: Gaussian -> Integer
norm (a :+ b) = (a ^ two) + (b ^ two)
  where
    two = 2 :: Integer
      -- This prevents compiler warnings about defaulting types.

conjugate :: Gaussian -> Gaussian
conjugate (a :+ b) = (a :+ negate b)
  -- The complex conjugate.
  -- Effectively flips about the real axis.

rotate_cw :: Gaussian -> Gaussian
rotate_cw (a :+ b) = (b :+ negate a)
  -- Rotate a gaussian integer clockwise.
  -- Equivalent to multiplying by i.

rotate_ccw :: Gaussian -> Gaussian
rotate_ccw (a :+ b) = (negate b :+ a)
  -- Rotate a gaussian integer counterclockwise.
  -- Equivalent to multiplying by -i.

flip_x :: Gaussian -> Gaussian
flip_x (a :+ b) = (negate a :+ b)
  -- Flips on imaginary axis.
  -- Analogous to conjugate.

flip_y :: Gaussian -> Gaussian
flip_y = conjugate
  -- Alias for conjugate.
  -- See conjugate.

swap_x_y :: Gaussian -> Gaussian
swap_x_y (a :+ b) = (b :+ a)
  -- Utility function, swaps real and imaginary axes.

toComplex :: Num p => Gaussian -> Complex.Complex p
toComplex (a :+ b) = (fromInteger a Complex.:+ fromInteger b)
  -- Converts a Gaussian to any complex number (Num p => Complex p)
  -- Fully qualified operators to avoid clashing.

As far as I know, my brace/indentation style isn't used anywhere else, so I'm a bit worried about that.

Answer 1

This looks good. All necessary type signatures are there, no type is too strict nor too generic. The indentation style is familiar, and the braces on data Gaussian are easy to read, so no worries about that.

However, there are a few oddities. Let's face them.

Num's signum

There aren't any laws for Num in the Haskell report, however there are some expectations. The operations on (+) and (*) should implement a ring and that's certainly true here. However, Num's documentation indicates that the following should hold:

abs x * signum x = x -- (1)

Unfortunately, we cannot take the easy way and just use

signum x = x / abs x

since the result is unlikely to be in Gaussian's domain. As we cannot hold (1), we're free to do whatever we want, but it should be reasonable. How does it become reasonable? By documenting your decision:

signum = id -- < your documentation here >

Prefix vs. postfix documentation

That took me by surprise. Most of Haskell's documentation is placed before the documented item, not after. The same holds for most other documentation in other languages, where the documentation is placed before or in the function itself:

/** @brief Calculates the absolute value of a Gaussian integer.
 *  @param x is a Gaussian integer
 *  @returns x's absolute value
*/
int g_abs(gint x) {
    return x.r * x.r + x.i * x.i;
}

def g_abs(x):
    """Returns the absolute value of x"""
    return x.r * x.r + x.i * x.i

The usual Haddock documentation also placed the documentation before the item, although it's possible to use postfix.

You're welcome to stay in the postfix style, but you should indicate that your documentation refers to the preceding item, for example with Haddock's ^:

i :: Gaussian
i = 0 :+ 1
--^ The imaginary unit, equal to sqrt(-1).

That removes possible ambiguity.

two = 2 :: Integer and x ^ two

The line two = 2 :: Integer violates the DRY principle, thus any optimization or bugfix on x ^ two needs to get repeated throughout the program, which is error prone.

However, there is a greater issue. GHC rewrites x ^ 2 to x * x (see Note [Inlining (^) in base's documention). An explicit variable might hinder that process.

Instead, you should provide a square function or just immediately multiply a value with itself. Both variants will get rid of the warnings and additional bindings, for example:

--| The Num instance for Gaussian.
--  We use the euclidean norm as abs. norm(a + bi) = a^2 + b^2
instance Num Gaussian
  where
    ...
    abs (a :+ b) = (a * a + b * b) :+ 0
    signum = id -- chosen at random
    ...

By the way: if you chose square x = x ^ (2 :: Int) you won't need to repeat your documentation.

Precision

You always want to stay as precise as possible, unless you are aiming for speed. Therefore, you should try to stay in exact types till the very last moment in magnitude.

Instead of

magnitude :: Floating c => Gaussian -> c
magnitude (a :+ b) = sqrt $ (square $ fromInteger a) + (square $ fromInteger b)

use

magnitude :: Floating c => Gaussian -> c
magnitude (a :+ b) = sqrt $ fromInteger (square a + square b)

You could also reuse abs if you split it into two functions:

instance Num Gaussian where
    ...
    abs x = (absIntegral x :+ 0)
    ...

absIntegral :: Gausian -> Integer
absIntegral (a :+ b) = square a + square b

magnitude :: Floating c => Gaussian -> c
magnitude = sqrt . fromInteger . absIntegral