Search code examples
haskellmatrixinverse

Inverse matrix in Haskell


I just started to learn Haskell.
And I have a question. I trying write a function to find the inverse matrix.
My type of matrix looks like that:

data Matrix a = M 
    { nrows :: !Int
    , ncols :: !Int
    , mvect :: V.Vector (V.Vector a)
    } deriving Eq

Also I have fromLists function.
The function for finding the determinant looks like that:

det :: Num a => Matrix a -> a
det (M 1 1 v) = V.head (V.head v)
det m =
    sum [ (-1)^(i-1) * m ! (i,1) * det (minorMatrix i 1 m) | i <- [1 .. nrows m] ]

So, my code for finding the inverse matrix:

coords :: Matrix a -> [[(Int, Int)]]
coords = zipWith (map . (,)) [0..] . map (zipWith const [0..])

delmatrix :: Int -> Int -> Matrix a -> Matrix a
delmatrix i j = dellist i . map (dellist j)
    where dellist i xs = take i xs ++ drop (i + 1) xs

mapMatrix :: (a -> a) -> Matrix a -> Matrix a
mapMatrix f = map (map f)

cofactor :: Num a => Int -> Int -> Matrix a -> a
cofactor i j m = ((-1) ** fromIntegral (i + j)) * det (delmatrix i j m)

cofactorM :: Matrix a -> Matrix a
cofactorM m = map (map (\(i,j) -> cofactor j i m)) $ coords m

inverse :: (Num a, Fractional a) => Matrix a -> Matrix a
inverse m = mapMatrix (* recip deter) $ cofactorM m
    where deter = det m

And what I have on a console:

Prelude> :r
[1 of 1] Compiling Matrixlab        ( Matrixlab.hs, interpreted )

Matrixlab.hs:120:38:
    Couldn't match expected type `Matrix a' with actual type `[a0]'
    Expected type: Matrix a -> [[Int]]
      Actual type: [a0] -> [b0]
    In the return type of a call of `map'
    In the second argument of `(.)', namely
      `map (zipWith const [0 .. ])'

Matrixlab.hs:123:17:
    Couldn't match expected type `Matrix a' with actual type `[a0]'
    Expected type: [a0] -> Matrix a
      Actual type: [a0] -> [a0]
    In the return type of a call of `dellist'
    In the first argument of `(.)', namely `dellist i'

Matrixlab.hs:127:15:
    Couldn't match expected type `Matrix a' with actual type `[a0]'
    Expected type: Matrix a -> Matrix a
      Actual type: [a0] -> [b0]
    In the return type of a call of `map'
    In the expression: map (map f)

Matrixlab.hs:130:24:
    Could not deduce (Floating a) arising from a use of `**'
    from the context (Num a)
      bound by the type signature for
                 cofactor :: Num a => Int -> Int -> Matrix a -> a
      at Matrixlab.hs:130:1-71
    Possible fix:
      add (Floating a) to the context of
        the type signature for
          cofactor :: Num a => Int -> Int -> Matrix a -> a
    In the first argument of `(*)', namely
      `((- 1) ** fromIntegral (i + j))'
    In the expression:
      ((- 1) ** fromIntegral (i + j)) * det (delmatrix i j m)
    In an equation for `cofactor':
        cofactor i j m
          = ((- 1) ** fromIntegral (i + j)) * det (delmatrix i j m)

Matrixlab.hs:133:15:
    Couldn't match expected type `Matrix a' with actual type `[[a]]'
    In the expression:
      map (map (\ (i, j) -> cofactor j i m)) $ coords m
    In an equation for `cofactorM':
        cofactorM m = map (map (\ (i, j) -> cofactor j i m)) $ coords m
Failed, modules loaded: none.

Help me, please.


Solution

  • First of all, if you want to use vectors, you should use vectors. If you are going to convert them back and from lists all the time, you should just stick with lists. The vector package provides all of the functions you have on lists for the Vector type. So you should be using those. In fact, it will almost always be easier to work with the Vector functions, since they give you built-in access to the index of an item. If you don't know what functions you have available, take a look at the docs.

    You also can't pretend that your Matrix is equivalent to a list of lists and expect the compiler to know what you're talking about. You have to take the vectors out of the Matrix constructor using pattern matching, operate on them, and put them back in. For example:

    import qualified Data.Vector as V
    
    cofactorM :: Floating a => Matrix a -> Matrix a
    cofactorM m@(M x y v) = M x y $ (V.imap $ \i -> V.imap $ \j -> const (cofactor i j m)) v 
    
    cofactor :: Floating a => Int -> Int -> Matrix a -> a
    cofactor ...
    

    Keep in mind that Vector is zero indexed, while your det function implies that your columns and rows are one indexed. You will have to compensate for this.

    Lastly, there is no good reason to package the number of rows and columns in with the Vector (Vector a) since getting the length of a vector is a constant time operation anyways.