(I added `Cc: ublas-dev' because I think that this can be interesting there, too.) Toon Knapen wrote:
On Wednesday 06 November 2002 16:41, Hickman, Greg wrote:
Why isn't there an overloaded
matrix_expression operator* (const matrix_expression&, const matrix_expression&);
for peforming matrix multiplication?
This issue has been discussed before and the trouble IIRC is that operator* could also mean element-wise matrix multiplication. IMHO the optimal solution would be to take over all conventions from Matlab, a convention most are familiar with, but the element-wise multiplication operator '.*' can't be used in C++ ;( Another suggestion was to defined the element-wise and the normal multiplication in a different namespace (I still like this idea very much)
Let's try: ====================================================== #include <iostream> #include <boost/numeric/ublas/matrix.hpp> #include <boost/numeric/ublas/io.hpp> namespace boost { namespace numeric { namespace ublas { namespace normal { template<class E1, class E2> BOOST_UBLAS_INLINE typename matrix_matrix_binary_traits< typename E1::value_type, E1, typename E2::value_type, E2 >::result_type operator * (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2) { return prod (e1, e2); } } }}} int main() { boost::numeric::ublas::matrix<double> m1 (2, 2), m2 (2, 2); // init m1 & m2 using namespace boost::numeric::ublas::normal; std::cout << m1 * m2 << std::endl; } ====================================================== g++ 3.2 complains (`slightly' edited ;o) : ====================================================== ambiguous overload for `matrix<double>& * matrix<double>&' operator candidates are: (defined here): matrix_matrix_binary_traits<>::result_type operator*(const matrix_expression<E>&, const matrix_expression<E2>&) [with E1 = matrix<double>, E2 = matrix<double>] (defined in `matrix_expression.hpp'): matrix_binary_scalar2_traits<>::result_type operator*(const matrix_expression<E>&, const T2&) [with E1 = matrix<double>, T2 = matrix<double>] (defined in `matrix_expression.hpp'): matrix_binary_scalar1_traits<>::result_type operator*(const T1&, matrix_expression<E2>&) [with T1 = matrix<double>, E2 = matrix<double>] ====================================================== Comeau's compiler 4.3.0.1 gives similar diagnostics. In `operator*(matrix_expression, T)' and `operator*(T, matrix_expression)' T is meant to be scalar, but compilers do not know this. Let's try to tell them (probably can be simplified): ==================================================== template<class E2> BOOST_UBLAS_INLINE typename matrix_binary_scalar1_traits< typename E2::value_type, E2, scalar_multiplies<typename E2::value_type, typename E2::value_type>
::result_type operator * (const typename E2::value_type &e1, const matrix_expression<E2> &e2) { // etc. ====================================================
It works; not only `m1 * m2', but also `m1 * trans (m2)', `m1 * (m1 + m2)' etc. Interesting thing is that if first parameter is declared as `const typename matrix_expression<E2>::value_type &e1' neither como nor g++ can find appropriate operator* for e.g. `2. * m1'. So, header `operator_mult.hpp' can be: ================================================= #include <boost/numeric/ublas/vector.hpp> #include <boost/numeric/ublas/matrix.hpp> namespace boost { namespace numeric { namespace ublas { namespace elemwise { // there are no elementwise multiplications in ublas, // at least not yet } namespace normal { template<class E1, class E2> BOOST_UBLAS_INLINE typename vector_scalar_binary_traits< E1, E2, vector_inner_prod< typename E1::value_type, typename E2::value_type, typename promote_traits< typename E1::value_type, typename E2::value_type >::promote_type > >::result_type operator * (const vector_expression<E1> &e1, const vector_expression<E2> &e2) { return inner_prod (e1, e2); } template<class E1, class E2> BOOST_UBLAS_INLINE typename matrix_vector_binary1_traits< typename E1::value_type, E1, typename E2::value_type, E2 >::result_type operator * (const matrix_expression<E1> &e1, const vector_expression<E2> &e2) { return prod (e1, e2); } template<class E1, class E2> BOOST_UBLAS_INLINE typename matrix_vector_binary2_traits< typename E1::value_type, E1, typename E2::value_type, E2 >::result_type operator * (const vector_expression<E1> &e1, const matrix_expression<E2> &e2) { return prod (e1, e2); } template<class E1, class E2> BOOST_UBLAS_INLINE typename matrix_matrix_binary_traits< typename E1::value_type, E1, typename E2::value_type, E2 >::result_type operator * (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2) { return prod (e1, e2); } } }}} ================================================ but, before it can be used, interfaces of existings operators * in matrix_expression.hpp and vector_expression.hpp must be adapted as described. There's one more, independant reason to change the interface of the operator* (scalar, Matrix_or_Vector_expression) and related/similar operators to something like operator * (const typename E2::value_type &e1, const M_or_V_expression<E2> &e2) Currently, one sometimes gets rather surprising and confusing results: ============================================== boost::numeric::ublas::matrix<double> m1 (2, 2); m1 (0, 0) = 0.1; m1 (0, 1) = 0.2; m1 (1, 0) = 0.3; m1 (1, 1) = 0.4; std::cout << 2 * m1 << std::endl; // note: `2', not `2.', std::cout << m1 * 2 << std::endl; // i.e. `int', not `double' ============================================== Output is: [2,2]((0,0),(0,0)) [2,2]((0.2,0.4),(0.6,0.8)) Sincerely, fres