ghc/interpreter/lib/Complex.hs

   1 -----------------------------------------------------------------------------
   2 -- Standard Library: Complex numbers
   3 --
   4 -- Suitable for use with Hugs 98
   5 -----------------------------------------------------------------------------
   6
   7 module Complex(Complex((:+)), realPart, imagPart, conjugate, mkPolar,
   8                cis, polar, magnitude, phase)  where
   9
  10 infix  6  :+
  11
  12 data (RealFloat a) => Complex a = !a :+ !a
  13                       deriving (Eq,Read,Show)
  14
  15 realPart, imagPart :: (RealFloat a) => Complex a -> a
  16 realPart (x:+y)     = x
  17 imagPart (x:+y)     = y
  18
  19 conjugate          :: (RealFloat a) => Complex a -> Complex a
  20 conjugate (x:+y)    = x :+ (-y)
  21
  22 mkPolar            :: (RealFloat a) => a -> a -> Complex a
  23 mkPolar r theta     = r * cos theta :+ r * sin theta
  24
  25 cis                :: (RealFloat a) => a -> Complex a
  26 cis theta           = cos theta :+ sin theta
  27
  28 polar              :: (RealFloat a) => Complex a -> (a,a)
  29 polar z             = (magnitude z, phase z)
  30
  31 magnitude, phase   :: (RealFloat a) => Complex a -> a
  32 magnitude (x:+y)    = scaleFloat k
  33                        (sqrt ((scaleFloat mk x)^2 + (scaleFloat mk y)^2))
  34                       where k  = max (exponent x) (exponent y)
  35                             mk = - k
  36 phase (0:+0)        = 0
  37 phase (x:+y)        = atan2 y x
  38
  39 instance (RealFloat a) => Num (Complex a) where
  40     (x:+y) + (x':+y')  = (x+x') :+ (y+y')
  41     (x:+y) - (x':+y')  = (x-x') :+ (y-y')
  42     (x:+y) * (x':+y')  = (x*x'-y*y') :+ (x*y'+y*x')
  43     negate (x:+y)      = negate x :+ negate y
  44     abs z              = magnitude z :+ 0
  45     signum 0           = 0
  46     signum z@(x:+y)    = x/r :+ y/r where r = magnitude z
  47     fromInteger n      = fromInteger n :+ 0
  48     fromInt n          = fromInt n :+ 0
  49
  50 instance (RealFloat a) => Fractional (Complex a) where
  51     (x:+y) / (x':+y')  = (x*x''+y*y'') / d :+ (y*x''-x*y'') / d
  52                          where x'' = scaleFloat k x'
  53                                y'' = scaleFloat k y'
  54                                k   = - max (exponent x') (exponent y')
  55                                d   = x'*x'' + y'*y''
  56     fromRational a     = fromRational a :+ 0
  57     fromDouble a       = fromDouble a :+ 0
  58
  59 instance (RealFloat a) => Floating (Complex a) where
  60     pi            = pi :+ 0
  61     exp (x:+y)    = expx * cos y :+ expx * sin y
  62                     where expx = exp x
  63     log z         = log (magnitude z) :+ phase z
  64     sqrt 0        = 0
  65     sqrt z@(x:+y) = u :+ (if y < 0 then -v else v)
  66                     where (u,v) = if x < 0 then (v',u') else (u',v')
  67                           v'    = abs y / (u'*2)
  68                           u'    = sqrt ((magnitude z + abs x) / 2)
  69     sin (x:+y)    = sin x * cosh y :+ cos x * sinh y
  70     cos (x:+y)    = cos x * cosh y :+ (- sin x * sinh y)
  71     tan (x:+y)    = (sinx*coshy:+cosx*sinhy)/(cosx*coshy:+(-sinx*sinhy))
  72                     where sinx  = sin x
  73                           cosx  = cos x
  74                           sinhy = sinh y
  75                           coshy = cosh y
  76     sinh (x:+y)   = cos y * sinh x :+ sin  y * cosh x
  77     cosh (x:+y)   = cos y * cosh x :+ sin y * sinh x
  78     tanh (x:+y)   = (cosy*sinhx:+siny*coshx)/(cosy*coshx:+siny*sinhx)
  79                     where siny  = sin y
  80                           cosy  = cos y
  81                           sinhx = sinh x
  82                           coshx = cosh x
  83     asin z@(x:+y) =  y':+(-x')
  84                      where  (x':+y') = log (((-y):+x) + sqrt (1 - z*z))
  85     acos z@(x:+y) =  y'':+(-x'')
  86                      where (x'':+y'') = log (z + ((-y'):+x'))
  87                            (x':+y')   = sqrt (1 - z*z)
  88     atan z@(x:+y) =  y':+(-x')
  89                      where (x':+y') = log (((1-y):+x) / sqrt (1+z*z))
  90     asinh z       = log (z + sqrt (1+z*z))
  91     acosh z       = log (z + (z+1) * sqrt ((z-1)/(z+1)))
  92     atanh z       = log ((1+z) / sqrt (1-z*z))
  93
  94 -----------------------------------------------------------------------------