{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE RebindableSyntax #-}

-- | The 'Normed' class represents elements of a normed vector
--   space. We define instances for all common numeric types.
module Normed
where

import BigFloat

import NumericPrelude hiding ( abs )
import Algebra.Absolute ( abs )
import qualified Algebra.Absolute as Absolute ( C )
import qualified Algebra.Algebraic as Algebraic ( C )
import Algebra.Algebraic ( root )
import qualified Algebra.RealField as RealField ( C )
import qualified Algebra.ToInteger as ToInteger ( C )
import qualified Algebra.ToRational as ToRational ( C )
import Data.Vector.Fixed ( S, Z )
import qualified Data.Vector.Fixed as V (
  Arity,
  map,
  maximum )
import Data.Vector.Fixed.Boxed ( Vec )

import Linear.Vector ( element_sum )


-- | Instances of the 'Normed' class know how to compute their own
--   p-norms for p=1,2,...,infinity.
--
class Normed a where
  norm_p :: (ToInteger.C c, Algebraic.C b, Absolute.C b) => c -> a -> b
  norm_infty :: (RealField.C b) => a -> b

  -- | The \"usual\" norm. Defaults to the 2-norm.
  norm :: (Algebraic.C b, Absolute.C b) => a -> b
  norm = norm_p (2 :: Integer)

-- Define instances for common numeric types.
instance Normed Integer where
  norm_p _ = abs . fromInteger
  norm_infty = abs . fromInteger

instance Normed Rational where
  norm_p _ = abs . fromRational'
  norm_infty = abs . fromRational'

instance Epsilon e => Normed (BigFloat e) where
  norm_p _ = abs . fromRational' . toRational
  norm_infty = abs . fromRational' . toRational

instance Normed Float where
  norm_p _ = abs . fromRational' . toRational
  norm_infty = abs . fromRational' . toRational

instance Normed Double where
  norm_p _ = abs . fromRational' . toRational
  norm_infty = abs . fromRational' . toRational


-- | 'Normed' instance for vectors of length zero. These are easy.
instance Normed (Vec Z a) where
  norm_p _ = const (fromInteger 0)
  norm_infty = const (fromInteger 0)


-- | 'Normed' instance for vectors of length greater than zero. We
--   need to know that the length is non-zero in order to invoke
--   V.maximum. We will generally be working with n-by-1 /matrices/
--   instead of vectors, but sometimes it's convenient to have these
--   instances anyway.
--
--   Examples:
--
--   >>> import Data.Vector.Fixed (mk3)
--   >>> import Linear.Vector (Vec3)
--   >>> let b = mk3 1 2 3 :: Vec3 Double
--   >>> norm_p 1 b :: Double
--   6.0
--   >>> norm b == sqrt 14
--   True
--   >>> norm_infty b :: Double
--   3.0
--
instance (V.Arity n, Absolute.C a, ToRational.C a, Ord a)
         => Normed (Vec (S n) a) where
  norm_p p x =
    (root p') $ element_sum $ V.map element_function x
    where
      element_function y = fromRational' $ (toRational y)^p'
      p' = toInteger p

  norm_infty x = fromRational' $ toRational $ V.maximum $ V.map abs x