From 3815bc25f62124063e02af83fe3c907336dc86f5 Mon Sep 17 00:00:00 2001
From: Alberto Ruiz <aruiz@um.es>
Date: Fri, 28 Sep 2007 12:25:56 +0000
Subject: algorithms interface reorganized

---
 HSSL.cabal                         |   1 +
 examples/tests.hs                  |   1 +
 lib/Data/Packed.hs                 |  70 +++++++++++++++++++++
 lib/Data/Packed/Internal/Matrix.hs |  12 ++++
 lib/Data/Packed/Internal/aux.h     |   2 +
 lib/Data/Packed/Matrix.hs          |   8 ++-
 lib/GSL/Matrix.hs                  |  35 +----------
 lib/GSL/gsl-aux.h                  |   2 -
 lib/LinearAlgebra.hs               |  16 +++--
 lib/LinearAlgebra/Algorithms.hs    | 125 ++++++++++++++++---------------------
 lib/LinearAlgebra/Interface.hs     |   7 +++
 lib/LinearAlgebra/Linear.hs        |  43 +++----------
 12 files changed, 171 insertions(+), 151 deletions(-)
 create mode 100644 lib/Data/Packed.hs

diff --git a/HSSL.cabal b/HSSL.cabal
index 22ac658..a3f4697 100644
--- a/HSSL.cabal
+++ b/HSSL.cabal
@@ -20,6 +20,7 @@ Exposed-modules:    Data.Packed.Internal,
                     Data.Packed.Internal.Common,
                     Data.Packed.Internal.Vector
                     Data.Packed.Internal.Matrix,
+                    Data.Packed,
                     Data.Packed.Vector,
                     Data.Packed.Matrix,
                     GSL.Vector,
diff --git a/examples/tests.hs b/examples/tests.hs
index b542113..5c14f96 100644
--- a/examples/tests.hs
+++ b/examples/tests.hs
@@ -328,6 +328,7 @@ tests = do
      [ test "arith1" $ ((ones (100,100) * 5 + 2)/0.5 - 7)**2 |~| (49 :: RM)
      , test "arith2" $ (((1+i) .* ones (100,100) * 5 + 2)/0.5 - 7)**2 |~| ( (140*i-51).*1 :: CM)
      , test "arith3" $ exp (i.*ones(10,10)*pi) + 1 |~| 0
+     , test "<\\>"   $ (3><2) [2,0,0,3,1,1::Double] <\> 3|>[4,9,5] |~| 2|>[2,3]
      ]
     putStrLn "--------- GSL ------"
     quickCheck $ \v -> ifft (fft v) |~| v
diff --git a/lib/Data/Packed.hs b/lib/Data/Packed.hs
new file mode 100644
index 0000000..668d2f7
--- /dev/null
+++ b/lib/Data/Packed.hs
@@ -0,0 +1,70 @@
+{-# OPTIONS_GHC -fglasgow-exts #-}
+-----------------------------------------------------------------------------
+{- |
+Module      :  Data.Packed
+Copyright   :  (c) Alberto Ruiz 2006-7
+License     :  GPL-style
+
+Maintainer  :  Alberto Ruiz (aruiz at um dot es)
+Stability   :  provisional
+Portability :  uses ffi
+
+The Vector and Matrix types and some utilities.
+
+-}
+-----------------------------------------------------------------------------
+
+module Data.Packed (
+    module Data.Packed.Vector,
+    module Data.Packed.Matrix,
+    module Data.Complex,
+    Container(..)
+) where
+
+import Data.Packed.Vector
+import Data.Packed.Matrix
+import Data.Complex
+import Data.Packed.Internal
+
+-- | conversion utilities
+class (Field e) => Container c e where
+    toComplex   :: RealFloat e => (c e, c e) -> c (Complex e)
+    fromComplex :: RealFloat e => c (Complex e) -> (c e, c e)
+    comp        :: RealFloat e => c e -> c (Complex e)
+    conj        :: RealFloat e => c (Complex e) -> c (Complex e)
+    real        :: c Double -> c e
+    complex     :: c e -> c (Complex Double)
+
+instance Container Vector Double where
+    toComplex = Data.Packed.Internal.toComplex
+    fromComplex = Data.Packed.Internal.fromComplex
+    comp = Data.Packed.Internal.comp
+    conj = Data.Packed.Internal.conj
+    real = id
+    complex = Data.Packed.comp
+
+instance Container Vector (Complex Double) where
+    toComplex = undefined -- can't match
+    fromComplex = undefined
+    comp = undefined
+    conj = undefined
+    real = Data.Packed.comp
+    complex = id
+
+instance Container Matrix Double where
+    toComplex = uncurry $ liftMatrix2 $ curry Data.Packed.toComplex
+    fromComplex z = (reshape c r, reshape c i)
+        where (r,i) = Data.Packed.fromComplex (cdat z)
+              c = cols z
+    comp = liftMatrix Data.Packed.Internal.comp
+    conj = liftMatrix Data.Packed.Internal.conj
+    real = id
+    complex = Data.Packed.comp
+
+instance Container Matrix (Complex Double) where
+    toComplex = undefined
+    fromComplex = undefined
+    comp = undefined
+    conj = undefined
+    real = Data.Packed.comp
+    complex = id
diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs
index 63ebddf..e76500b 100644
--- a/lib/Data/Packed/Internal/Matrix.hs
+++ b/lib/Data/Packed/Internal/Matrix.hs
@@ -24,6 +24,8 @@ import Complex
 import Control.Monad(when)
 import Data.List(transpose,intersperse)
 import Data.Maybe(fromJust)
+import Foreign.C.String
+import Foreign.C.Types
 
 -----------------------------------------------------------------
 
@@ -371,6 +373,16 @@ fromComplex z = (r,i) where
 comp :: Vector Double -> Vector (Complex Double)
 comp v = toComplex (v,constant 0 (dim v))
 
+-- | loads a matrix efficiently from formatted ASCII text file (the number of rows and columns must be known in advance).
+fromFile :: FilePath -> (Int,Int) -> IO (Matrix Double)
+fromFile filename (r,c) = do
+    charname <- newCString filename
+    res <- createMatrix RowMajor r c
+    c_gslReadMatrix charname // mat dat res // check "gslReadMatrix" []
+    --free charname  -- TO DO: free the auxiliary CString
+    return res
+foreign import ccall "aux.h matrix_fscanf" c_gslReadMatrix:: Ptr CChar -> TM
+
 -------------------------------------------------------------------------
 
 -- Generic definitions
diff --git a/lib/Data/Packed/Internal/aux.h b/lib/Data/Packed/Internal/aux.h
index 83111e5..73334e3 100644
--- a/lib/Data/Packed/Internal/aux.h
+++ b/lib/Data/Packed/Internal/aux.h
@@ -26,3 +26,5 @@ int diagR(KRVEC(d),RMAT(r));
 int diagC(KCVEC(d),CMAT(r));
 
 const char * gsl_strerror (const int gsl_errno);
+
+int matrix_fscanf(char*filename, RMAT(a));
diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs
index 404fde7..2b93348 100644
--- a/lib/Data/Packed/Matrix.hs
+++ b/lib/Data/Packed/Matrix.hs
@@ -25,9 +25,10 @@ module Data.Packed.Matrix (
     fromBlocks,
     flipud, fliprl,
     subMatrix, takeRows, dropRows, takeColumns, dropColumns,
+    extractRows,
     ident, diag, diagRect, takeDiag,
     liftMatrix, liftMatrix2,
-    format, dispR, readMatrix, fromArray2D
+    format, dispR, readMatrix, fromFile, fromArray2D
 ) where
 
 import Data.Packed.Internal
@@ -209,3 +210,8 @@ dispC d m = disp m (shfc d)
 -- | creates a matrix from a table of numbers.
 readMatrix :: String -> Matrix Double
 readMatrix = fromLists . map (map read). map words . filter (not.null) . lines
+
+-- | rearranges the rows of a matrix according to the order given in a list of integers. 
+extractRows :: Field t => [Int] -> Matrix t -> Matrix t
+extractRows l m = fromRows $ extract (toRows $ m) l
+    where extract l is = [l!!i |i<-is]
diff --git a/lib/GSL/Matrix.hs b/lib/GSL/Matrix.hs
index c1bce37..3a54226 100644
--- a/lib/GSL/Matrix.hs
+++ b/lib/GSL/Matrix.hs
@@ -18,17 +18,14 @@ module GSL.Matrix(
     qr,
     cholR, -- cholC,
     luSolveR, luSolveC,
-    luR, luC,
-    fromFile, extractRows
+    luR, luC
 ) where
 
 import Data.Packed.Internal
 import Data.Packed.Matrix(fromLists,ident,takeDiag)
 import GSL.Vector
 import Foreign
-import Foreign.C.Types
 import Complex
-import Foreign.C.String
 
 {- | eigendecomposition of a real symmetric matrix using /gsl_eigen_symmv/.
 
@@ -257,7 +254,7 @@ p is a permutation:
 
 L \* U obtains a permuted version of the original matrix:
 
-@\> 'extractRows' p m
+@\> extractRows p m
   2.+3.i   -7.   0.
       1.    2.  -3.
       1.  -1.i  2.i
@@ -298,35 +295,7 @@ luC m = (l,u,p, fromIntegral s') where
     add = liftMatrix2 $ vectorZipC Add
     mul = liftMatrix2 $ vectorZipC Mul
 
-extract l is = [l!!i |i<-is]
-
 {- auxiliary function to get triangular matrices
 -}
 triang r c h v = reshape c $ fromList [el i j | i<-[0..r-1], j<-[0..c-1]]
     where el i j = if j-i>=h then v else 1 - v
-
-{- | rearranges the rows of a matrix according to the order given in a list of integers. 
-
-> > extractRows [3,3,0,1] (ident 4)
-> 0. 0. 0. 1.
-> 0. 0. 0. 1.
-> 1. 0. 0. 0.
-> 0. 1. 0. 0.
-
--}
-extractRows :: Field t => [Int] -> Matrix t -> Matrix t
-extractRows l m = fromRows $ extract (toRows $ m) l
-
---------------------------------------------------------------
-
--- | loads a matrix efficiently from formatted ASCII text file (the number of rows and columns must be known in advance).
-fromFile :: FilePath -> (Int,Int) -> IO (Matrix Double)
-fromFile filename (r,c) = do
-    charname <- newCString filename
-    res <- createMatrix RowMajor r c
-    c_gslReadMatrix charname // mat dat res // check "gslReadMatrix" []
-    --free charname  -- TO DO: free the auxiliary CString
-    return res
-foreign import ccall "gsl-aux.h matrix_fscanf" c_gslReadMatrix:: Ptr CChar -> TM
-
----------------------------------------------------------------------------
diff --git a/lib/GSL/gsl-aux.h b/lib/GSL/gsl-aux.h
index dd1a247..3ccac25 100644
--- a/lib/GSL/gsl-aux.h
+++ b/lib/GSL/gsl-aux.h
@@ -52,8 +52,6 @@ int integrate_qags(double f(double,void*), double a, double b, double prec, int
 
 int polySolve(KRVEC(a), CVEC(z));
 
-int matrix_fscanf(char*filename, RMAT(a));
-
 int minimize(double f(int, double*), double tolsize, int maxit, 
                  KRVEC(xi), KRVEC(sz), RMAT(sol));
 
diff --git a/lib/LinearAlgebra.hs b/lib/LinearAlgebra.hs
index d8f44b4..d6ef647 100644
--- a/lib/LinearAlgebra.hs
+++ b/lib/LinearAlgebra.hs
@@ -13,17 +13,15 @@ Basic matrix computations implemented by BLAS, LAPACK and GSL.
 -}
 -----------------------------------------------------------------------------
 module LinearAlgebra (
-    module Data.Packed.Vector,
-    module Data.Packed.Matrix,
-    module LinearAlgebra.Instances,
-    module LinearAlgebra.Interface,
+    module Data.Packed,
+    module LinearAlgebra.Linear,
     module LinearAlgebra.Algorithms,
-    module Complex
+    module LinearAlgebra.Instances,
+    module LinearAlgebra.Interface
 ) where
 
+import Data.Packed
+import LinearAlgebra.Linear
+import LinearAlgebra.Algorithms
 import LinearAlgebra.Instances
 import LinearAlgebra.Interface
-import LinearAlgebra.Algorithms
-import Data.Packed.Matrix
-import Data.Packed.Vector
-import Complex
\ No newline at end of file
diff --git a/lib/LinearAlgebra/Algorithms.hs b/lib/LinearAlgebra/Algorithms.hs
index 162426f..a67f822 100644
--- a/lib/LinearAlgebra/Algorithms.hs
+++ b/lib/LinearAlgebra/Algorithms.hs
@@ -9,106 +9,114 @@ Maintainer  :  Alberto Ruiz (aruiz at um dot es)
 Stability   :  provisional
 Portability :  uses ffi
 
-A simple interface to the available matrix computations. We have defined some generic functions 
-on both real and complex matrices in such a way that higher level algorithms and 
-testing properties are valid for both base types.
+A generic interface for a number of essential functions. Using it some higher level algorithms
+and testing properties can be written for both real and complex matrices.
 
-In any case the specific functions for particular base types can also be explicitly
+In any case, the specific functions for particular base types can also be explicitly
 imported from the LAPACK and GSL.Matrix modules.
 
 -}
 -----------------------------------------------------------------------------
 
 module LinearAlgebra.Algorithms (
--- * Basic Linear Algebra
-    scale,
-    add,
-    multiply, dot, outer,
-    linearSolve, linearSolveSVD,
+-- * Linear Systems
+    linearSolve,
+    inv, pinv,
+    pinvTol, det,
 -- * Matrix factorizations
-    svd, lu, eig, eigSH,
-    qr, chol,
--- * Utilities
-    Normed(..), NormType(..),
-    det,inv,pinv,full,economy,
-    pinvTol,
---    pinvTolg,
+-- ** Singular value decomposition
+    svd,
+    full, economy,
+-- ** Eigensystems
+    eig, LinearAlgebra.Algorithms.eigS, LinearAlgebra.Algorithms.eigH,
+-- ** Other 
+    LinearAlgebra.Algorithms.qr, chol,
+-- * Nullspace
     nullspacePrec,
     nullVector,
--- * Conversions
-    toComplex, fromComplex,
-    comp, real, complex,
-    conj, ctrans,
 -- * Misc
     eps, i,
-    scaleRecip,
-    addConstant,
-    sub,
-    mul,
-    divide,
+    ctrans,
+    Normed(..), NormType(..),
+    GenMat(linearSolveSVD,lu,eigSH)
 ) where
 
 
 import Data.Packed.Internal hiding (fromComplex, toComplex, comp, conj)
-import Data.Packed.Matrix
+import Data.Packed
 import GSL.Matrix
 import GSL.Vector
 import LAPACK
 import Complex
 import LinearAlgebra.Linear
 
--- | Base types for which some optimized matrix computations are available
-class (Linear Matrix t) => Optimized t where
+-- | matrix computations available for both real and complex matrices
+class (Linear Matrix t) => GenMat t where
     svd         :: Matrix t -> (Matrix t, Vector Double, Matrix t)
     lu          :: Matrix t -> (Matrix t, Matrix t, [Int], t)
     linearSolve :: Matrix t -> Matrix t -> Matrix t
     linearSolveSVD :: Matrix t -> Matrix t -> Matrix t
-    ctrans :: Matrix t -> Matrix t
     eig         :: Matrix t -> (Vector (Complex Double), Matrix (Complex Double))
     eigSH       :: Matrix t -> (Vector Double, Matrix t)
     chol        :: Matrix t -> Matrix t
+    -- | conjugate transpose
+    ctrans :: Matrix t -> Matrix t
 
-instance Optimized Double where
+instance GenMat Double where
     svd = svdR
     lu  = luR
     linearSolve = linearSolveR
     linearSolveSVD = linearSolveSVDR Nothing
     ctrans = trans
     eig = eigR
-    eigSH = eigS
+    eigSH = LAPACK.eigS
     chol = cholR
 
-instance Optimized (Complex Double) where
+instance GenMat (Complex Double) where
     svd = svdC
     lu  = luC
     linearSolve = linearSolveC
     linearSolveSVD = linearSolveSVDC Nothing
     ctrans = conjTrans
     eig = eigC
-    eigSH = eigH
+    eigSH =  LAPACK.eigH
     chol = error "cholC not yet implemented" -- waiting for GSL-1.10
 
+-- | eigensystem of a symmetric matrix
+eigS :: Matrix Double -> (Vector Double, Matrix Double)
+eigS = LAPACK.eigS
+
+-- | eigensystem of a hermitian matrix
+eigH :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
+eigH = LAPACK.eigH
+
+qr :: Matrix Double -> (Matrix Double, Matrix Double)
+qr = GSL.Matrix.qr
+
 square m = rows m == cols m
 
-det :: Optimized t => Matrix t -> t
+det :: GenMat t => Matrix t -> t
 det m | square m = s * (product $ toList $ takeDiag $ u)
       | otherwise = error "det of nonsquare matrix"
     where (_,u,_,s) = lu m
 
-inv :: Optimized t => Matrix t -> Matrix t
+inv :: GenMat t => Matrix t -> Matrix t
 inv m | square m = m `linearSolve` ident (rows m)
       | otherwise = error "inv of nonsquare matrix"
 
-pinv :: Optimized t => Matrix t -> Matrix t
+pinv :: GenMat t => Matrix t -> Matrix t
 pinv m = linearSolveSVD m (ident (rows m))
 
-
+full :: Field t 
+     => (Matrix t -> (Matrix t, Vector Double, Matrix t)) -> Matrix t -> (Matrix t, Matrix Double, Matrix t)
 full svd m = (u, d ,v) where
     (u,s,v) = svd m
     d = diagRect s r c
     r = rows m
     c = cols m
 
+economy :: Field t 
+        => (Matrix t -> (Matrix t, Vector Double, Matrix t)) -> Matrix t -> (Matrix t, Vector Double, Matrix t)
 economy svd m = (u', subVector 0 d s, v') where
     (u,s,v) = svd m
     sl@(g:_) = toList (complex s)
@@ -123,39 +131,25 @@ economy svd m = (u', subVector 0 d s, v') where
     v' = takeColumns d v
 
 
-{- | Machine precision of a Double.
-
->> eps
-> 2.22044604925031e-16
-
-(The value used by GNU-Octave)
-
--}
+-- | The machine precision of a Double: @eps == 2.22044604925031e-16@ (the value used by GNU-Octave).
 eps :: Double
 eps =  2.22044604925031e-16
 
-{- | The imaginary unit
-
-@> 'ident' 3 \<\> i
-1.i   0.   0.
- 0.  1.i   0.
- 0.   0.  1.i@
-
--}
+-- | The imaginary unit: @i == 0.0 :+ 1.0@
 i :: Complex Double
 i = 0:+1
 
 
 -- | matrix product
-mXm :: (Num t, Optimized t) => Matrix t -> Matrix t -> Matrix t
+mXm :: (Num t, GenMat t) => Matrix t -> Matrix t -> Matrix t
 mXm = multiply
 
 -- | matrix - vector product
-mXv :: (Num t, Optimized t) => Matrix t -> Vector t -> Vector t
+mXv :: (Num t, GenMat t) => Matrix t -> Vector t -> Vector t
 mXv m v = flatten $ m `mXm` (asColumn v)
 
 -- | vector - matrix product
-vXm :: (Num t, Optimized t) => Vector t -> Matrix t -> Vector t
+vXm :: (Num t, GenMat t) => Vector t -> Matrix t -> Vector t
 vXm v m = flatten $ (asRow v) `mXm` m
 
 
@@ -167,15 +161,6 @@ norm2 = toScalarR Norm2
 norm1 :: Vector Double -> Double
 norm1 = toScalarR AbsSum
 
-vectorMax :: Vector Double -> Double
-vectorMax = toScalarR Max
-vectorMin :: Vector Double -> Double
-vectorMin = toScalarR Min
-vectorMaxIndex :: Vector Double -> Int
-vectorMaxIndex = round . toScalarR MaxIdx
-vectorMinIndex :: Vector Double -> Int
-vectorMinIndex = round . toScalarR MinIdx
-
 data NormType = Infinity | PNorm1 | PNorm2 -- PNorm Int
 
 pnormRV PNorm2 = norm2
@@ -199,7 +184,7 @@ pnormCM Infinity m = vectorMax $ liftMatrix (liftVector magnitude) m `mXv` const
 --pnormCM _ _ = error "p norm not yet defined"
 
 -- -- | computes the p-norm of a matrix or vector (with the same definitions as GNU-octave). pnorm 0 denotes \\inf-norm. See also 'norm'.
---pnorm :: (Container t, Optimized a) => Int -> t a -> Double
+--pnorm :: (Container t, GenMat a) => Int -> t a -> Double
 --pnorm = pnormG
 
 class Normed t where
@@ -222,7 +207,7 @@ instance Normed (Matrix (Complex Double)) where
 -----------------------------------------------------------------------
 
 -- | The nullspace of a matrix from its SVD decomposition.
-nullspacePrec :: Optimized t
+nullspacePrec :: GenMat t
               => Double          -- ^ relative tolerance in 'eps' units
               -> Matrix t   -- ^ input matrix
               -> [Vector t] -- ^ list of unitary vectors spanning the nullspace
@@ -235,12 +220,12 @@ nullspacePrec t m = ns where
     ns = drop rank $ toRows $ ctrans v
 
 -- | The nullspace of a matrix, assumed to be one-dimensional, with default tolerance (shortcut for @last . nullspacePrec 1@).
-nullVector :: Optimized t => Matrix t -> Vector t
+nullVector :: GenMat t => Matrix t -> Vector t
 nullVector = last . nullspacePrec 1
 
 ------------------------------------------------------------------------
 
-{- | Pseudoinverse of a real matrix with the desired tolerance, expressed as a
+{-  Pseudoinverse of a real matrix with the desired tolerance, expressed as a
 multiplicative factor of the default tolerance used by GNU-Octave (see 'pinv').
 
 @\> let m = 'fromLists' [[1,0,    0]
@@ -258,7 +243,7 @@ multiplicative factor of the default tolerance used by GNU-Octave (see 'pinv').
 0. 0. 1.@
 
 -}
-pinvTol :: Double -> Matrix Double -> Matrix Double
+--pinvTol :: Double -> Matrix Double -> Matrix Double
 pinvTol t m = v' `mXm` diag s' `mXm` trans u' where
     (u,s,v) = svdR m
     sl@(g:_) = toList s
diff --git a/lib/LinearAlgebra/Interface.hs b/lib/LinearAlgebra/Interface.hs
index 3392d54..0c65a8b 100644
--- a/lib/LinearAlgebra/Interface.hs
+++ b/lib/LinearAlgebra/Interface.hs
@@ -16,6 +16,7 @@ Operators for frequent operations.
 
 module LinearAlgebra.Interface(
     (<>),(<.>),
+    (<\>),
     (.*),(*/),
     (<|>),(<->),
 ) where
@@ -23,6 +24,7 @@ module LinearAlgebra.Interface(
 import LinearAlgebra.Linear
 import Data.Packed.Vector
 import Data.Packed.Matrix
+import LinearAlgebra.Algorithms
 
 class Mul a b c | a b -> c where
  infixl 7 <>
@@ -59,6 +61,11 @@ a .* x = scale a x
 infixl 7 */
 v */ x = scale (recip x) v
 
+-- | least squares solution of a linear system, similar to the \\ operator of Matlab\/Octave (based on linearSolveSVD).
+(<\>) :: (GenMat a) => Matrix a -> Vector a -> Vector a
+infixl 7 <\>
+m <\> v = flatten (linearSolveSVD m (reshape 1 v))
+
 ------------------------------------------------
 
 class Joinable a b where
diff --git a/lib/LinearAlgebra/Linear.hs b/lib/LinearAlgebra/Linear.hs
index 85daa4a..3e6c55d 100644
--- a/lib/LinearAlgebra/Linear.hs
+++ b/lib/LinearAlgebra/Linear.hs
@@ -9,10 +9,10 @@ Maintainer  :  Alberto Ruiz (aruiz at um dot es)
 Stability   :  provisional
 Portability :  uses ffi
 
+Basic optimized operations on vectors and matrices.
 
 -}
 -----------------------------------------------------------------------------
--- #hide
 
 module LinearAlgebra.Linear (
     Linear(..),
@@ -21,25 +21,22 @@ module LinearAlgebra.Linear (
 
 
 import Data.Packed.Internal
-import Data.Packed.Matrix
+import Data.Packed
 import GSL.Vector
 import Complex
 
 -- | basic optimized operations
-class (Field e) => Linear c e where
+class (Container c e) => Linear c e where
     scale       :: e -> c e -> c e
-    scaleRecip  :: e -> c e -> c e
     addConstant :: e -> c e -> c e
     add         :: c e -> c e -> c e
     sub         :: c e -> c e -> c e
+    -- | element by element multiplication
     mul         :: c e -> c e -> c e
+    -- | element by element division
     divide      :: c e -> c e -> c e
-    toComplex   :: RealFloat e => (c e, c e) -> c (Complex e)
-    fromComplex :: RealFloat e => c (Complex e) -> (c e, c e)
-    comp        :: RealFloat e => c e -> c (Complex e)
-    conj        :: RealFloat e => c (Complex e) -> c (Complex e)
-    real        :: c Double -> c e
-    complex     :: c e -> c (Complex Double)
+    -- | scale the element by element reciprocal of the object: @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
+    scaleRecip  :: e -> c e -> c e
 
 instance Linear Vector Double where
     scale = vectorMapValR Scale
@@ -49,12 +46,6 @@ instance Linear Vector Double where
     sub = vectorZipR Sub
     mul = vectorZipR Mul
     divide = vectorZipR Div
-    toComplex = Data.Packed.Internal.toComplex
-    fromComplex = Data.Packed.Internal.fromComplex
-    comp = Data.Packed.Internal.comp
-    conj = Data.Packed.Internal.conj
-    real = id
-    complex = LinearAlgebra.Linear.comp
 
 instance Linear Vector (Complex Double) where
     scale = vectorMapValC Scale
@@ -64,12 +55,6 @@ instance Linear Vector (Complex Double) where
     sub = vectorZipC Sub
     mul = vectorZipC Mul
     divide = vectorZipC Div
-    toComplex = undefined -- can't match
-    fromComplex = undefined
-    comp = undefined
-    conj = undefined
-    real = LinearAlgebra.Linear.comp
-    complex = id
 
 instance Linear Matrix Double where
     scale x = liftMatrix (scale x)
@@ -79,14 +64,6 @@ instance Linear Matrix Double where
     sub = liftMatrix2 sub
     mul = liftMatrix2 mul
     divide = liftMatrix2 divide
-    toComplex = uncurry $ liftMatrix2 $ curry LinearAlgebra.Linear.toComplex
-    fromComplex z = (reshape c r, reshape c i)
-        where (r,i) = LinearAlgebra.Linear.fromComplex (cdat z)
-              c = cols z
-    comp = liftMatrix Data.Packed.Internal.comp
-    conj = liftMatrix Data.Packed.Internal.conj
-    real = id
-    complex = LinearAlgebra.Linear.comp
 
 instance Linear Matrix (Complex Double) where
     scale x = liftMatrix (scale x)
@@ -96,12 +73,6 @@ instance Linear Matrix (Complex Double) where
     sub = liftMatrix2 sub
     mul = liftMatrix2 mul
     divide = liftMatrix2 divide
-    toComplex = undefined
-    fromComplex = undefined
-    comp = undefined
-    conj = undefined
-    real = LinearAlgebra.Linear.comp
-    complex = id
 
 --------------------------------------------------
 
-- 
cgit v1.2.3