Re: [boost] [rfc] A light weight statically typed matrix template class (again?)

From: Chris Fairles (chris.fairles_at_[hidden])
Date: 2007-12-06 14:48:37

You might be interested in looking at blitz++ or tvmet. They both use
expression templates to statically unroll various operators (for small
sized vectors/matricies anyway). If I ever find a spare minute, I've
always wanted to boostify one of those libs.

Boost also has a UBLAS impl. which has a matrix and vector class
(although not static).

Cheers,
Chris

On Dec 6, 2007 1:37 PM, Andreas Harnack <ah.boost.04_at_[hidden]> wrote:
> Hi @,
>
> shouldn't there be at least an elementary matrix class available
> in the standard? I notice there have been some discussions in
> the past, but there doesn't seem to be anything in the pipeline.
> So, is there any interest in a (really!) light weight statically
> typed matrix template class? I'm thinking of something like:
>
>
> #include <algorithm>
> #include <functional>
>
> template<typename T, unsigned int M, unsigned int N>
> class matrix
> {
>      protected:
>          T m[M][N];
>      public:
>          typedef T value_type;
>          enum { rows = M, columns = N };
>          T* data() { return &**m; }
>          T const* data() const { return &**m; }
>
>          matrix() {}
>          matrix(T const& c) {
>              std::fill(data(), data()+M*N, c); }
>          matrix(matrix const& m) {
>              std::copy(m.data(), m.data()+M*N, data()); }
>
>          template<typename U>
>          matrix(matrix<U,M,N> const& m) {
>              std::copy(m.data(), m.data()+M*N, data()); }
>
>          template<typename U>
>          matrix& operator=(matrix<U,M,N> const& m) {
>              std::copy(m.data(),m.data()+M*N,data());
>              return *this; }
>
>          template<class UnaryFunction>
>          matrix& map(UnaryFunction);
>          template<class BinaryFunction>
>          matrix& join(matrix const&, BinaryFunction);
>
>          matrix& operator+=(matrix const& m) {
>              return join(m, std::plus<T>() ); }
>          matrix& operator-=(matrix const& m) {
>              return join(m, std::minus<T>() ); }
>          matrix& operator*=(T const& c) {
>              return map(std::bind2nd(std::multiplies<T>(), c)); }
>          matrix& operator/=(T const& c) {
>              return map(std::bind2nd(std::divides<T>(), c)); }
>          matrix& operator%=(T const& c) {
>              return map(std::bind2nd(std::modulus<T>(), c)); }
>
>          T& operator()(unsigned int i, unsigned int j)
>              { return m[i][j]; }
>          T const& operator()(unsigned int i, unsigned int j)
>              const { return m[i][j]; }
> };
>
>
> template<typename T, unsigned int M, unsigned int N>
>      template<class UnaryFunction> matrix<T,M,N>&
>      matrix<T,M,N>::map( UnaryFunction op) {
>          std::transform(data(), data()+M*N, data(), op);
>          return *this; }
>
> template<typename T, unsigned int M, unsigned int N>
>      template<class BinaryFunction> matrix<T,M,N>&
>      matrix<T,M,N>::join(matrix const& m, BinaryFunction op) {
>          std::transform(data(), data()+M*N, m.data(),
>              data(), op); return *this; }
>
> // global converter; allows type deduction
> template<typename T, unsigned int M, unsigned int N> inline
>      matrix<T,M,N> mtrx(T (&a)[M][N]) { matrix<T,M,N> tmp;
>          std::copy(&**a, &**a+M*N, tmp.data()); return tmp; }
>
> // unary and binary operators
> template<typename T, unsigned int M, unsigned int N> inline
>      bool operator==(matrix<T,M,N> const& m1, matrix<T,M,N>
>          const& m2) { return std::equal(
>              m1.data(), m1.data()+M*N, m2.data()); }
>
> template<typename T, unsigned int M, unsigned int N> inline
>      bool operator!=(matrix<T,M,N> const& m1, matrix<T,M,N>
>          const& m2) { return ! std::equal(
>              m1.data(), m1.data()+M*N, m2.data()); }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline
>      operator-(matrix<T,M,N> const& m) {
>          matrix<T,M,N> tmp = T(); return tmp-=m; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator+(matrix<T,M,N> tmp,
>      matrix<T,M,N> const& m) { return tmp+=m; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator-(matrix<T,M,N> tmp,
>      matrix<T,M,N> const& m) { return tmp-=m; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator*(matrix<T,M,N> tmp,
>      T const& c) { return tmp*=c; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator*(T const& c,
>      matrix<T,M,N> tmp) { return tmp*=c; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator/(matrix<T,M,N> tmp,
>      T const& c) { return tmp/=c; }
>
> template<typename T, unsigned int M, unsigned int N>
> matrix<T,M,N> inline operator%(matrix<T,M,N> tmp,
>      T const& c) { return tmp%=c; }
>
>
> The advantage of this approach is its simplicity, in fact, it's
> so straight forward I'm surprised it's not already in the
> standard. There's nothing in there (yet), that's not already in
> the STL. It's just like bundling a few things together and
> making them available in a more convenient notation  and/or from
> a different point of view.
>
> There are still a few essential but fairly basic things missing,
> like matrix multiplication and transposing a matrix. Elementary
> row operations might be useful as well for those, who want to
> provide complex algorithms, but that's it, basically. There
> shouldn't be anything in there that's more complicated than
> that. (Matrices are such a general concept that it can be very
> tricky, if not impossible, to provide any guaranties for
> anything above the basic stuff.)
>
> Such a matrix class would be useful in all situations, where
> matrix notation can lead to shorter and better readable code,
> while matrices remain reasonably small and/or performance is not
> a primary concern. I'm thinking of domains like graphical
> programming. Matrices could be passed around by value and almost
> be treated as build in types. They could also be used as
> building block for more advanced libraries.
>
> I came across the problem of a missing matrix class when
> refining my Euclid class. There are some application problems,
> that are best solved using matrices (like vector space
> transformation). The current release therefore contains a
> (diagonal) matrix class, but I feel it doesn't really belong
> there. Matrices are by far general enough, to live in there own
> library. The available Linear Algebra libraries, however, are by
> far too complex for such simple tasks.
>
> I've uploaded a first version in the vault (simple_matrix.zip),
> to have something to start with. It already contains some
> missing features in a very naive implementation and comes with a
> test suit and an example program.
>
> Any interest? Comments most welcome!
> Andreas
>
>
> _______________________________________________
> Unsubscribe & other changes: http://lists.boost.org/mailman/listinfo.cgi/boost
>