Dear OCamlers,

In 1994, Barton and Nackman in their book 'Scientific Engineering in
C++' [1] demonstrated how one could encode the rules of Dimensional
Analysis [2] into the C++ type system enabling compile-time checking
(no runtime-cost) of the plausibility (at least up to the dimensional
correctness) of computations.

In 2004, Abrahams & Gurtovy in 'C++ Template Metaprogramming' [3]
showed the Barton Nackman technique to be elegantly implementable
using compile time type sequences encoding integer constants. At the
end of this post, I provide a complete listing of their example
program [4].

The key properties of the system (as I see it) are:
  - Encoding of integers as types; 
  - Compile time manipulation of sequences of these integer encodings
    to deduce/produce new derived types.

Now, it is not immediately obvious to me how to approach this problem
in OCaml. It irks me some that I can't immediately produce a yet more
elegant OCaml program for this problem and leaves me feeling like C++
has "got something over on us" here ;)

My question therefore is: Does anyone have suggestions/pointers
on how to approach automatic dimensional analysis via the OCaml type
system? 

Best,

-- 
Shayne Fletcher

[1] John J. Barton and Lee R. Nackman. Scientific and Engineering C++:
    an Introduction with Advanced Techniques and Examples. Reading,
    MA: Addison Wesley. ISBN 0-201-53393-6. 1994.


[3] David Abrahams and Aleksey Gurtovy C++ Template Metaprogramming:
    Concepts, Tools, and Techniques from Boost and Beyond (C++ in
    Depth Series), Addison-Wesley Professional. ISBN:0321227255. 2004.

[4] Code listing:

    //"c:/program files (x86)/Microsoft Visual Studio 10.0/vc/vcvarsall.bat" x64
    //cl /Fedimension.exe /EHsc /I d:/boost_1_55_0 dimension.cpp
    
    #include <boost/mpl/vector_c.hpp>
    #include <boost/mpl/transform.hpp>
    #include <boost/mpl/placeholders.hpp>
    #include <boost/mpl/equal.hpp>
    #include <boost/mpl/plus.hpp>
    #include <boost/mpl/minus.hpp>
    #include <boost/static_assert.hpp>
    
    typedef boost::mpl::vector_c<int,1,0,0,0,0,0,0> mass;
    typedef boost::mpl::vector_c<int,0,1,0,0,0,0,0> length;
    typedef boost::mpl::vector_c<int,0,0,1,0,0,0,0> time;
    typedef boost::mpl::vector_c<int,0,0,0,1,0,0,0> charge;
    typedef boost::mpl::vector_c<int,0,0,0,0,1,0,0> temperature;
    typedef boost::mpl::vector_c<int,0,0,0,0,0,1,0> intensity;
    typedef boost::mpl::vector_c<int,0,0,0,0,0,0,1> angle;
    typedef boost::mpl::vector_c<int,0,1,-1,0,0,0,0> velocity;     // l/t
    typedef boost::mpl::vector_c<int,0,1,-2,0,0,0,0> acceleration; // l/(t2)
    typedef boost::mpl::vector_c<int,1,1,-1,0,0,0,0> momentum;     // ml/t
    typedef boost::mpl::vector_c<int,1,1,-2,0,0,0,0> force;        // ml/(t2)
    typedef boost::mpl::vector_c<int,0,0,0,0,0,0,0> scalar;
    
    template <class T, class Dimensions>
    class quantity
    {
    public:
      explicit quantity (T val) 
        : val (val)
      {}
      template <class OtherDimensions>
      quantity (quantity<T, OtherDimensions> const& other)
        : val (other.value ()) {
        BOOST_MPL_ASSERT( (boost::mpl::equal<Dimensions, OtherDimensions>));
      }
      T value () const { return val; }
    private:
      T val;
    };
    
    template <class T, class D>
    quantity<T, D>
    operator + (quantity<T, D> x, quantity<T, D> y )
    {
      return quantity<T, D>(x.value () + y.value ());
    }
    
    template <class T, class D>
    quantity<T, D>
    operator - (quantity<T, D> x, quantity<T, D> y )
    {
      return quantity<T, D>(x.value () - y.value ());
    }
    
    template <class T, class D1, class D2>
    quantity <
      T
    , typename boost::mpl::transform<
        D1, D2, boost::mpl::plus<
                  boost::mpl::placeholders::_1
                , boost::mpl::placeholders::_2> >::type 
    >
    operator* (quantity<T, D1> x, quantity <T, D2> y)
    {
      typedef typename boost::mpl::transform<
        D1, D2, boost::mpl::plus<
                  boost::mpl::placeholders::_1
                  , boost::mpl::placeholders::_2> >::type D;
    
      return quantity<T, D> (x.value () * y.value ());
    }
    
    template <class T, class D1, class D2>
    quantity <
      T
    , typename boost::mpl::transform<
        D1, D2, boost::mpl::minus<
                  boost::mpl::placeholders::_1
                , boost::mpl::placeholders::_2> >::type 
    >
    operator/ (quantity<T, D1> x, quantity <T, D2> y)
    {
      typedef typename boost::mpl::transform<
        D1, D2, boost::mpl::minus<
                  boost::mpl::placeholders::_1
                  , boost::mpl::placeholders::_2> >::type D;
    
      return quantity<T, D> (x.value () / y.value ());
    }
    
    // -- test
    
    #include <iostream>
    #include <limits>
    #include <cassert>
    
    int main ()
    {
      quantity<float, mass> m (5.0f);
      quantity<float, acceleration> a(9.8f);
      quantity<float, force> f = m * a;
      quantity<float, mass> m2 = f / a;
    
      assert ((std::abs ((m2 - m).value ())) <= std::numeric_limits<double>::epsilon ());
    
      return 0;
    }