Alex Rivera | Logout

Can good type systems distinguish between matrices in different bases?

Asked 2011-05-01T18:24:09.200
30

My program (Hartree-Fock/iterative SCF) has two matrices F and F' which are really the same matrix expressed in two different bases. I just lost three hours of debugging time because I accidentally used F' instead of F. In C++, the type-checker doesn't catch this kind of error because both variables are Eigen::Matrix<double, 2, 2> objects.

I was wondering, for the Haskell/ML/etc. people, whether if you were writing this program you would have constructed a type system where F and F' had different types? What would that look like? I'm basically trying to get an idea how I can outsource some logic errors onto the type checker.

Edit: The basis of a matrix is like the unit. You can say 1L or however many gallons, they both mean the same thing. Or, to give a vector example, you can say (0,1) in Cartesian coordinates or (1,pi/2) in polar. But even though the meaning is the same, the numerical values are different.

Edit: Maybe units was the wrong analogy. I'm not looking for some kind of record type where I can specify that the first field will be litres and the second gallons, but rather a way to say that this matrix as a whole, is defined in terms of some other matrix (the basis), where the basis could be any matrix of the same dimensions. E.g., the constructor would look something like mkMatrix [[1, 2], [3, 4]] [[5, 6], [7, 8]] and then adding that object to another matrix would type-check only if both objects had the same matrix as their second parameters. Does that make sense?

Edit: definition on Wikipedia, worked examples

haskell types ocaml sml type-systems
Edit
Report

1 Answer

20

This is entirely possible in Haskell.

Statically checked dimensions

Haskell has arrays with statically checked dimensions, where the dimensions can be manipulated and checked statically, preventing indexing into the wrong dimension. Some examples:

This will only work on 2-D arrays:

multiplyMM :: Array DIM2 Double -> Array DIM2 Double -> Array DIM2 Double

An example from repa should give you a sense. Here, taking a diagonal requires a 2D array, returns a 1D array of the same type.

diagonal :: Array DIM2 e -> Array DIM1 e

or, from Matt sottile's repa tutorial, statically checked dimensions on a 3D matrix transform:

f :: Array DIM3 Double -> Array DIM2 Double
f u =
  let slabX = (Z:.All:.All:.(0::Int))
      slabY = (Z:.All:.All:.(1::Int))
      u' = (slice u slabX) * (slice u slabX) +
           (slice u slabY) * (slice u slabY)
  in
    R.map sqrt u'

Statically checked units

Another example from outside of matrix programming: statically checked units of dimension, making it a type error to confuse e.g. feet and meters, without doing the conversion.

 Prelude> 3 *~ foot + 1 *~ metre
 1.9144 m

or for a whole suite of SI units and quanitie

answered 2011-05-01T18:36:09.243

Your Answer