how do I pass the boost matrix prod() function as a multiplies function?

Question

I'm trying to perform a matrix exponentiation, but I don't want to copy/paste my exponentiation function, and would rather use class templates. The problem is that for boost matrices, to multiply matrices, you use the prod function (instead of operator*).

It seems that g++ isn't able to figure out the template I want to use. The error I'm getting with the below code is

41:37: error: no matching function for call to 'my_pow(boost::numeric::ublas::matrix<int>&, int, <unresolved overloaded function type>)'

here's the code:

#include <iostream>
using namespace std;

#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>

typedef long long int64;

template <class T, class M> T my_pow(T b, int64 e, M mult)
{
  if (e == 1) return b;
  if (e % 2 == 1) return mult(b, my_pow(b, e - 1, mult));
  T tmp = my_pow(b, e / 2, mult);
  return mult(tmp, tmp);
}
template <class T> T my_pow(T b, int64 e) { return my_pow(b, e, multiplies<T>()); }

int main()
{
  using namespace boost::numeric::ublas;
  matrix<int> m(3, 3);
  for (unsigned i = 0; i < m.size1(); ++i)
    for (unsigned j = 0; j < m.size2(); ++j)
      m(i, j) = 3 * i + j;
  std::cout << m << std::endl;
  std::cout << my_pow(m, 2, prod) << std::endl;
}

Is there any way to pass prod() to my_pow so the template resolves? Thanks.

It case it's not clear: b is the base, e is the exponent, and my_pow is to compute b^e

Answer 1

The reason why you are getting the compiler error is that there are many overloads of the prod function and at the call to my_pow the compiler needs to know which one to provide. The compiler cannot deduce that you will be applying the pow function to the first argument of your function, so it is at a loss here.

One solution would be to explicitly cast the function pointer to the right type, but for the uBlas prod overloads determining the right type to cast to can be quite complex.

Another solution is to create a polymorphic function object that delegates to the appropriate pow function. Note that the implementation below makes the huge assumption that prod( m, m) returns a value of the same type as m (or something convertible to it), but then again, this is the same assumption that your my_pow makes and the temporaries that this creates are hard to avoid if the power e can only be determined at run-time.

An example of a polymorphic function class that would do the trick:

struct my_prod
{
    template< typename M>
    M operator()( const M &left, const M &right) const
    {
        return prod( left, right);
    }

};

Now, if you change your call to my_pow into this:

std::cout << my_pow(m, 2, my_prod()) << std::endl;

It should work (it does for me).

Answer 2

There are two problems. First, prod is a templatized function, so you can't just pass prod as a function pointer. Instead you would need to pass prod<...> with the specific template paramaters filled in.

However in this case that still won't fix your issue because even with the specified template parameters, prod still has several overloads and the compiler cannot determine which one it should use. It is possible to fix this by declaring a function pointer that specifies the arguments and return type. However since ublas uses complex template metaprogramming this will be extremely ugly and I would not recommend it. Instead I would write a wrapper function around prod to call the specific overload you want. Here is a pretty generic wrapper that should work with any ublas matrix:

template <class E1, class E2> 
typename boost::numeric::ublas::matrix_matrix_binary_traits<
        typename E1::value_type, E1, 
        typename E2::value_type, E2>::result_type
my_prod(const boost::numeric::ublas::matrix_expression<E1>& e1, 
        const boost::numeric::ublas::matrix_expression<E1>& e2)
{
    return prod(e1, e2);
}

Then you can call my_pow using my_prod like so:

my_pow(m, 2, my_prod<matrix<int>, matrix<int> >)

---

Just for fun though, here is the function pointer declaration you would need to pass to resolve the template paramaters and overloads. This declares a function pointer named prod_ptr that points to the the specific overload of prod you want:

matrix_matrix_binary_traits<matrix<int>::value_type, matrix<int>, matrix<int>::value_type, matrix<int> >::result_type 
    (*prod_ptr)(const matrix_expression<matrix<int> >&, const matrix_expression<matrix<int> >&) = 
    &prod<matrix_matrix_binary_traits<matrix<int>::value_type, matrix<int>, matrix<int>::value_type, matrix<int> >::result_type, matrix<int>, matrix<int> >;

Then you would be able to call my_pow using the function pointer:

my_pow(m, 2, prod_ptr);