{-# LANGUAGE DeriveDataTypeable, GeneralizedNewtypeDeriving, StandaloneDeriving, FlexibleInstances, FlexibleContexts, MultiParamTypeClasses, Rank2Types, UndecidableInstances #-}
{- |
Module      :  Numeric.VariablePrecision.Complex
Copyright   :  (c) Claude Heiland-Allen 2012
License     :  BSD3

Maintainer  :  claude@mathr.co.uk
Stability   :  unstable
Portability :  DeriveDataTypeable, GeneralizedNewtypeDeriving, StandaloneDeriving, FlexibleInstances, FlexibleContexts, MultiParamTypeClasses, Rank2Types, UndecidableInstances

Complex numbers with variable precision.

-}
module Numeric.VariablePrecision.Complex
  ( module Data.Complex.Generic
  , VComplex()
  , toComplex
  , fromComplex
  , withComplex
  , recodeComplex
  , scaleComplex
  , scaleVComplex
  , toComplexDFloat
  , toComplexDFixed
  , fromComplexDFloat
  , fromComplexDFixed
  , withComplexDFloat
  , withComplexDFixed
  ) where

import Data.Complex.Generic
import Data.Complex.Generic.Default

import Numeric.VariablePrecision.Fixed
import Numeric.VariablePrecision.Float
import Numeric.VariablePrecision.Precision
import Numeric.VariablePrecision.Algorithms (recodeFloat)

-- | Newtype wrapper around 'Complex' so that instances can be written
--   for 'HasPrecision' and 'VariablePrecision'.
newtype VComplex t p = FromComplex
  { -- | Convert 'VComplex' to 'Complex'.
    toComplex :: Complex (t p)
  }
  deriving (Eq)

deriving instance Num (Complex (t p)) => Num (VComplex t p)
deriving instance Fractional (Complex (t p)) => Fractional (VComplex t p)
deriving instance Floating (Complex (t p)) => Floating (VComplex t p)

-- | Convert 'Complex' to 'VComplex'.
fromComplex :: Complex (t p) -> VComplex t p
fromComplex = FromComplex

instance NaturalNumber p => Num (Complex (VFloat p)) where
  (+) = addDefault
  (-) = subDefault
  (*) = mulDefault
  negate = negateDefault
  fromInteger = fromIntegerDefault
  abs = absDefault
  signum = signumDefault

instance NaturalNumber p => Num (Complex (VFixed p)) where
  (+) = addDefault
  (-) = subDefault
  (*) = mulDefault
  negate = negateDefault
  fromInteger = fromIntegerDefault
  abs = error "Num.abs: unimplemented for Complex VFixed"
  signum = error "Num.signum: unimplemented for Complex VFixed"

instance NaturalNumber p => Fractional (Complex (VFloat p)) where
  (/) = divDefaultRF
  fromRational = fromRationalDefault

instance NaturalNumber p => Fractional (Complex (VFixed p)) where
  (/) = divDefault
  fromRational = fromRationalDefault

instance NaturalNumber p => Floating (Complex (VFloat p)) where
  pi = piDefault
  exp = expDefault
  log = logDefault
  sqrt = sqrtDefault
  sin = sinDefault
  cos = cosDefault
  tan = tanDefault
  sinh = sinhDefault
  cosh = coshDefault
  tanh = tanhDefault
  asin = asinDefault
  acos = acosDefault
  atan = atanDefault
  asinh = asinhDefault
  acosh = acoshDefault
  atanh = atanhDefault

instance NaturalNumber p => ComplexRect (Complex (VFloat p)) (VFloat p) where
  mkRect x y = (x :+ y)
  rect (x :+ y) = (x, y)
  real = realDefault
  imag = imagDefault
  realPart = realPartDefault
  imagPart = imagPartDefault
  conjugate = conjugateDefault
  magnitudeSquared = magnitudeSquaredDefault
  sqr = sqrDefaultRF
  (.*) = rmulDefault
  (*.) = mulrDefault

instance NaturalNumber p => ComplexRect (Complex (VFixed p)) (VFixed p) where
  mkRect x y = (x :+ y)
  rect (x :+ y) = (x, y)
  real = realDefault
  imag = imagDefault
  realPart = realPartDefault
  imagPart = imagPartDefault
  conjugate = conjugateDefault
  magnitudeSquared = magnitudeSquaredDefault
  sqr = sqrDefault -- FIXME
  (.*) = rmulDefault
  (*.) = mulrDefault

instance NaturalNumber p => ComplexPolar (Complex (VFloat p)) (VFloat p) where
  mkPolar = mkPolarDefault
  cis = cisDefault
  polar = polarDefault
  magnitude = magnitudeDefaultRF
  phase = phaseDefaultRF

instance NaturalNumber p => ComplexRect (VComplex VFloat p) (VFloat p) where
  mkRect x y = FromComplex (x :+ y)
  rect = rect . toComplex
  real = realDefault
  imag = imagDefault
  realPart = realPartDefault
  imagPart = imagPartDefault
  conjugate = conjugateDefault
  magnitudeSquared = magnitudeSquaredDefault
  sqr = sqrDefaultRF
  (.*) = rmulDefault
  (*.) = mulrDefault

instance NaturalNumber p => ComplexRect (VComplex VFixed p) (VFixed p) where
  mkRect x y = FromComplex (x :+ y)
  rect = rect . toComplex
  real = realDefault
  imag = imagDefault
  realPart = realPartDefault
  imagPart = imagPartDefault
  conjugate = conjugateDefault
  magnitudeSquared = magnitudeSquaredDefault
  sqr = sqrDefault -- FIXME
  (.*) = rmulDefault
  (*.) = mulrDefault

instance NaturalNumber p => ComplexPolar (VComplex VFloat p) (VFloat p) where
  mkPolar = mkPolarDefault
  cis = cisDefault
  polar = polarDefault
  magnitude = magnitudeDefaultRF
  phase = phaseDefaultRF

instance HasPrecision (VComplex VFloat)

instance HasPrecision (VComplex VFixed)

instance VariablePrecision (VComplex VFloat) where
  adjustPrecision = withComplex (fmap adjustPrecision)

instance VariablePrecision (VComplex VFixed) where
  adjustPrecision = withComplex (fmap adjustPrecision)

instance Normed (VComplex VFloat) where
  norm1 z = abs (realPart z) + abs (imagPart z)
  norm2 = magnitude
  norm2Squared = magnitudeSquared
  normInfinity z = abs (realPart z) `max` abs (imagPart z)

instance NaturalNumber p => Show (VComplex VFloat p) where
  showsPrec p = showsPrec p . toComplex

instance NaturalNumber p => Show (VComplex VFixed p) where
  showsPrec p = showsPrec p . toComplex

instance NaturalNumber p => Read (VComplex VFloat p) where
  readsPrec p = map (first fromComplex) . readsPrec p
    where first f (a, b) = (f a, b)

instance NaturalNumber p => Read (VComplex VFixed p) where
  readsPrec p = map (first fromComplex) . readsPrec p
    where first f (a, b) = (f a, b)

-- | Lift an operation on 'X.Complex' to one on 'VComplex'.
withComplex :: (Complex (t p) -> Complex (t q)) -> (VComplex t p -> VComplex t q)
withComplex f = fromComplex . f . toComplex

-- | Much like 'mapComplex' 'recodeFloat'.
recodeComplex :: (RealFloat a, RealFloat b) => Complex a -> Complex b
recodeComplex = fmap recodeFloat

-- | Much like 'withComplex' 'scaleComplex''.
scaleVComplex :: NaturalNumber p => Int -> VComplex VFloat p -> VComplex VFloat p
scaleVComplex = withComplex . scaleComplex

-- | Much like 'mapComplex' 'scaleFloat'.
scaleComplex :: RealFloat r => Int -> Complex r -> Complex r
scaleComplex = fmap . scaleFloat

-- | Freeze a 'VComplex VFloat'.
toComplexDFloat :: NaturalNumber p => VComplex VFloat p -> Complex DFloat
toComplexDFloat = fmap toDFloat . toComplex

-- | Freeze a 'VComplex VFixed'.
toComplexDFixed :: NaturalNumber p => VComplex VFixed p -> Complex DFixed
toComplexDFixed = fmap toDFixed . toComplex

-- | Thaw a 'Complex DFloat'.  Results in 'Nothing' on precision mismatch.
fromComplexDFloat :: NaturalNumber p => Complex DFloat -> Maybe (VComplex VFloat p)
fromComplexDFloat (x :+ y) = do
  r <- fromDFloat x
  i <- fromDFloat y
  return (r .+ i)

-- | Thaw a 'Complex DFixed'.  Results in 'Nothing' on precision mismatch.
fromComplexDFixed :: NaturalNumber p => Complex DFixed -> Maybe (VComplex VFixed p)
fromComplexDFixed (x :+ y) = do
  r <- fromDFixed x
  i <- fromDFixed y
  return (r .+ i)

-- | Thaw a 'Complex DFloat' to its natural precision.  'Nothing' is passed on
--   precision mismatch between real and imaginary parts.
withComplexDFloat :: Complex DFloat -> (forall p . NaturalNumber p => Maybe (VComplex VFloat p) -> r) -> r
withComplexDFloat (x :+ y) f = withDFloat x $ \r -> f $ do
  i <- fromDFloat y
  return (r .+ i)

-- | Thaw a 'Complex DFixed' to its natural precision.  'Nothing' is passed on
--   precision mismatch between real and imaginary parts.
withComplexDFixed :: Complex DFixed -> (forall p . NaturalNumber p => Maybe (VComplex VFixed p) -> r) -> r
withComplexDFixed (x :+ y) f = withDFixed x $ \r -> f $ do
  i <- fromDFixed y
  return (r .+ i)