C++

Forum

 
Operator overloading with matrices

cmajor28 (48)
So I have been creating a matrix class and I'm almost done with it. All I have left are the + - * operators (there is no matrix division). I already have += -= *= but I can't seem to get the syntax for the others. I'm using templates.
I want to know the syntax for a member function and a non-member function. I've been reading for an hour on the internet and I keep trying what I see and it doesn't work.
Duthomhas (13206)
1
2
3
4
5
 
template <...stuff...>
matrix <...stuff...> operator + ( const matrix <...stuff...> & a, const matrix <...stuff...> & b )
{
   ...
}

Hope this helps.
cmajor28 (48)
 
 
template <class U> Matrix<T> operator + (const Matrix<T> &lmat, const Matrix<U> &rmat);


Then what's wrong with this? (Note my class template is <class T>). I also tried: 

 
 
template <class U, class V> Matrix<U> operator + (const Matrix<U> &lmat, const Matrix<V> &rmat);


so that it wouldn't conflict with my class template.

EDIT: Okay so apparently I HAVE TO put the word "friend." So it looks like this now:
 
 
	template <class U> friend Matrix<T> operator + (const Matrix<T> &lmat, const Matrix<U> &rmat);
Last edited on
JLBorges (13770)
> I already have += -= *=

This is the canonical C++ way:

a. overload operators +, - and * as non-member non-friend functions
b. which take the first parameter by value
c. and delegate to the appropriate compound-assignment member function 


1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
 
template < typename T > struct matrix
{
    // ...

    matrix<T>& operator += ( const matrix<T>& that ) ;
    matrix<T>& operator -= ( const matrix<T>& that ) ;
    matrix<T>& operator *= ( const matrix<T>& that ) ;

    // ...
};

// note: the first parameter is *passed by value*
// note: this allows us to implement the binary arithmetic operator in terms of the compound assignment operator (less code to maintain)
// note: this also allows optimisation of chained operations like a+b+c-d

template < typename T > matrix<T> operator+ ( matrix<T> a, const matrix<T>& b ) 
{ 
    return a += b ; // delegate to compound assignment; return by value
}

template < typename T > matrix<T> operator- ( matrix<T> a, const matrix<T>& b ) { return a -= b ; }

template < typename T > matrix<T> operator* ( matrix<T> a, const matrix<T>& b ) { return a *= b ; }
cmajor28 (48)
Ahh! I see that makes so much more sense. I was wondering why I thought i needed friend. Thanks a lot!

EDIT: Well I only seem to be getting this to work if I do cout << A+B; If I do C = A + B I get a stack dump error.
Last edited on
JLBorges (13770)
> If I do C = A + B I get a stack dump error.

That seems to be a problem with your overloaded simple assignment operator. operator=
cmajor28 (48)
Well here are the two operators I was using.

Assignment
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
 
template <class T>
template <class U> Matrix<T>& operator = (const Matrix<U> &mat) {
	
	if (this != reinterpret_cast<Matrix<T>*>(const_cast<Matrix<U>*>(&mat))) {
		// Free memory
		free(matrix);
			
		// Set new rows and columns
		m = mat.row();
		n = mat.col();
		
		// Allocate memory for matrix
		matrix = (T*)malloc(m*n*sizeof(T));
		
		if (sizeof(T) != sizeof(U)) {
			size_t i,j;
			U* point = mat.data();
			
			// Set matrix to other matrix
			for (i=0;i<m;i++) {
				for (j=0;j<n;j++) {
					*(matrix+n*i+j) = *(point+n*i+j);
				}
			}
		}
		else {
			// Set matrix to copy of mat
			memcpy(matrix, mat.data(), m*n*sizeof(T));
		}
	}
	return *this;
}



Addition Assignment
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 
template <class T>
template <class U> Matrix<T>& Matrix<T>::operator += (const Matrix<U> &mat) {
	
	if (m != mat.row() || n != mat.col()) {
		clear();
	}
	else {
		size_t i,j;
		U* temp = mat.data();

		// Add corresponding elements
		for (i=0;i<m;i++) {
			for (j=0;j<n;j++) {
				*(matrix+n*i+j) += *(temp+n*i+j);
			}
		}
	}
	return *this;
}
JLBorges (13770)
Use std::vector<> perhaps?
And explicitly default the foundation operations - copy constructor, move constructor, copy assignment and move assignment.

Something along these lines, may be?

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
 
#include <iostream>
#include <vector>
#include <initializer_list>
#include <type_traits>
#include <algorithm>
#include <iomanip>

template < typename T > struct matrix
{
    matrix() = default ;

    matrix( std::size_t nrows, std::size_t ncols, const T& v = T() ) : data( nrows, std::vector<T>( ncols, v ) ) {}

    template < typename U > matrix( std::initializer_list< std::initializer_list<U> > ilist )
    {
        std::size_t max_size = 0 ;
        for( const auto& row : ilist )
        {
            data.emplace_back(row) ;
            max_size = std::max( max_size, data.back().size() ) ;
        }

        for( auto& row : data ) row.resize(max_size) ;
    }

    matrix( const matrix<T>& that ) = default ; // defaulted copy constructor
    matrix( matrix<T>&& that ) noexcept = default ; // defaulted move constructor

    template < typename U > matrix( const matrix<U>& that,
                                    const typename std::enable_if< std::is_convertible<U,T>::value >::type* = nullptr )
    {
        const auto& other_data = that.pub_data_helper() ;
        data.resize( other_data.size() ) ;
        for( std::size_t i = 0 ; i < other_data.size() ; ++i )
        {
            data[i].resize( other_data[i].size() ) ;
            std::copy( std::begin(other_data[i]), std::end(other_data[i]), std::begin(data[i]) ) ;
        }
    }

    matrix<T>& operator= ( const matrix<T>& that ) = default ; // defaulted copy assignment
    matrix<T>& operator= ( matrix<T>&& that ) noexcept = default ; // defaulted move assignment

    template < typename U >
    typename std::enable_if< std::is_convertible<U,T>::value, matrix<T>& >::type  operator= ( const matrix<U>& that )
    { return *this = matrix<T>(that) ; } // not the most efficient

    // implicitly declared destructor

    private: std::vector< std::vector<T> > data ;
    public: const std::vector< std::vector<T> >& pub_data_helper() const { return data ; }
};

template < typename C, typename TRAITS, typename T >
std::basic_ostream<C,TRAITS>& operator<< ( std::basic_ostream<C,TRAITS>& stm, const matrix<T>& m )
{
   const auto& data = m.pub_data_helper() ;
   for( const auto& row : data )
   {
       for( const T& v : row ) stm << v << stm.widen(' ') ;
       stm << stm.widen('\n') ;
   }
   return stm ;
}

int main()
{
    matrix<int> m1 { { {0,1,2,3}, {4,5,6,7} } };

    matrix<int> m2(m1) ;
    std::cout << m2 << '\n' ;

    matrix<double> m3(m1) ;
    std::cout << std::fixed << std::setprecision(1) << m3 << '\n' ;

    m2 = matrix<int> { { {0,1,2}, {3,4,5}, {6,7,8} } };
    std::cout << m2 << '\n' ;

    m3 = m2 ;
    std::cout << m3 << '\n' ;
}

http://coliru.stacked-crooked.com/a/e7d549cdd4cec5a3
cmajor28 (48)
I don't really want to use std::vector because it causes my constructor speed to be worse as well as how I perform my calculations. I initially was going to use vectors for everything but I did some speed test and using memory allocation with pointers was a little faster in my constructors and A LOT faster in computing determinants and inverses of matrices.

EDIT: well 1d vectors were good for a constructor but 2d vectors were bad! But you did give me a great idea! I like how you had support for initializer lists. I know I want to implement that. As for my problem I'm going to take it one step at a time until I find the issue.
Last edited on
Duthomhas (13206)
I don't see anything wrong with your current functions... but I could be missing something.

And you're right -- vector is slow.

You could, however, use some <algorithm>s to speed things up over a loop-loop:

15
16
17
 
		if (sizeof(T) != sizeof(U)) {
			std::copy( mat.data(), mat.data() + m*n, data() );
		}

Hope this helps.
cmajor28 (48)
Thats true! I might do that.
ne555 (10692)
1
2
 
//sizeof(T) == sizeof(U)
memcpy(matrix, mat.data(), m*n*sizeof(T));
in my system sizeof(int) = sizeof(float) = 4


> well 1d vectors were good for a constructor but 2d vectors were bad!
then use 1d vectors std::vector<T> matrix( rows*columns );
JLBorges (13770)
cmajor28> well 1d vectors were good for a constructor but 2d vectors were bad!
ne555>> then use 1d vectors std::vector<T> matrix( rows*columns );

ne555 +1


std::vector<> is not slower than dynamically allocated c-style arrays.
(Many years ago, vector implementations used to be measurably slower).

http://www.cplusplus.com/forum/general/126092/#msg683239
ln -s /Archive2/fd/2449d5921648c7/main.cpp foo.cpp 
clang++ -std=c++11 -stdlib=libc++ -O3 -Wall -pedantic-errors -c foo.cpp
clang++ -std=c++11 -stdlib=libc++ -O3 -Wall -pedantic-errors main.cpp foo.o -lsupc++
./a.out
g++-4.8 -std=c++11 -O3 -Wall -pedantic-errors -c foo.cpp
g++-4.8 -std=c++11 -O3 -Wall -pedantic-errors main.cpp foo.o
./a.out
vector: 1200 millisecs.
carray: 1190 millisecs.
vector: 1180 millisecs.
carray: 1200 millisecs.


http://www.cplusplus.com/forum/general/119436/#msg652470
g++-4.8 -std=c++11 -O3 -march=core2 -Wall -pedantic-errors main.cpp && ./a.out
vector: 3.06 secs.
c style array: 3.01 secs.


If the library has a good std::valarray<> implementation, consider using it.
Compound assignment operators, slices and friends make it a very attractive choice.
http://www.cplusplus.com/forum/general/151083/#msg787620
clang++ -std=c++11 -stdlib=libc++ -O3 -Wall -Wextra -pedantic-errors main.cpp -lsupc++ && ./a.out
echo --------- && g++ -std=c++11 -O3 -Wall -Wextra -pedantic-errors main.cpp && ./a.out
array: 130 msecs.
vector: 100 msecs.
valarray: 130 msecs.
1.26754e+08
---------
array: 140 msecs.
vector: 90 msecs.
valarray: 100 msecs.
4.85773e+07
Last edited on
cmajor28 (48)
Well when I tried computing the inverse of a 10x10 matrix, a carray was a good bit faster. With the c array, I was using pointer arithmetic which is just slightly faster than indexing but not enough to matter. However, in the matrix determinant algorithm (which I used the recursive expansion of minors because it works in every case and is pretty fast) the vector was a good bit slower because expansions of minors you must re-size the matrix through each recursion and using pointer arithmetic I could get by without doing that. But ne555 may have a great idea. If I can't get this to work I will definitely use that idea. I think I could make the expansion of minors algorithm to be just as fast with a 1d vector. And thanks for telling me that size of int == size of float. Any ideas on a workaround for what I was doing? Or does the memcpy function automatically change it?

EDIT: I actually found a way. if (typeid(int)==typeid(float)) returns false so I'll use that.

As for my error, my addition operator is perfect. It is definitely my simple assignment operator as JLBorges said. For some reason it does not work AT ALL. I cannot even say A = B;

Also note: somebody said something about a copy and move assignment operator. This is the only assignment operator I have for two matrices so should I implement the other type?

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
 
template <class T>
template <class U> Matrix<T>& Matrix<T>::operator = (const Matrix<U> &mat) {

	if (this != reinterpret_cast<Matrix<T>*>(const_cast<Matrix<U>*>(&mat))) {
		// Free memory
		free(matrix);
			
		// Set new rows and columns
		m = mat.row();
		n = mat.col();
		
		// Allocate memory for matrix
		matrix = (T*)malloc(m*n*sizeof(T));
		
		if (typeid(T)!=typeid(U)) {
			size_t i,j;
			U* point = mat.data();
			
			// Set matrix to other matrix
			for (i=0;i<m;i++) {
				for (j=0;j<n;j++) {
					*(matrix+n*i+j) = *(point+n*i+j);
				}
			}
		}
		else {
			// Set matrix to copy of mat
			memcpy(matrix, mat.data(), m*n*sizeof(T));
		}
	}
	return *this;
}


FINAL EDIT: I decided to separate it into 2 functions. One for the same data type and one for different. I like it a lot better to.
Last edited on
JLBorges (13770)
> I decided to separate it into 2 functions. One for the same data type and one for different. I like it a lot better too.

Yes, yes.


> Also note: somebody said something about a copy and move assignment operator.

Standard library components like vector are Moveable and MoveAssignable.
We don't need to do anything special to support move semantics if we use a vector.


> the vector was a good bit slower because expansions of minors you must re-size the matrix
> through each recursion and using pointer arithmetic I could get by without doing that.

std::vector<>::iterator is a random access iterator; it supports the complete range of pointer semantics.

It is trivial to access the elements of a std::vector<> via c-style raw pointers.
1
2
 
std::vector<int> vec = { 0, 1, 2, 3, 4, 5, 6, 7 } ;
const int* c_style_pointer_to_first_element = vec.empty() ? nullptr : &vec.front() ;



> But ne555 may have a great idea. If I can't get this to work I will definitely use that idea.
> I think I could make the expansion of minors algorithm to be just as fast with a 1d vector.

Yes. Consider using that idea even if you can get this to work.

The earlier code, modified to use a flattened std::vector<T> as suggested by ne555:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
 
#include <iostream>
#include <vector>
#include <initializer_list>
#include <type_traits>
#include <algorithm>
#include <iomanip>
#include <stdexcept>

template < typename T > struct matrix
{
    using iterator = typename std::vector<T>::iterator ;
    using const_iterator = typename std::vector<T>::const_iterator ;

    std::size_t rows() const { return nrows ; }
    std::size_t cols() const { return ncols ; }
    bool empty() const { return rows() == 0 ; }

    const std::vector<T>& data() const { return data_ ; } // low level access to stored elements

    // even lower level access to stored elements as a c-style array
    const T* c_ptr() const { return empty() ? nullptr : std::addressof(data_.front()) ; }

    struct row
    {
        iterator begin() const { return begin_ ; }
        iterator end() const { return end_ ; }
        std::size_t size() const { return end() - begin() ; }

        T& operator[] ( std::size_t pos ) { return begin()[pos] ; }
        const T& operator[] ( std::size_t pos ) const { return begin()[pos] ; }

        T& at( std::size_t pos )
        { if( pos >= size() ) throw std::out_of_range("matrix<T>::row: out of range element access") ; return begin()[pos] ; }
        const T& at( std::size_t pos ) const
        { if( pos >= size() ) throw std::out_of_range("matrix<T>::row: out of range element access") ; return begin()[pos] ; }

        const iterator begin_ ;
        const iterator end_ ;
    };

    struct const_row
    {
        const_iterator begin() const { return begin_ ; }
        const_iterator end() const { return end_ ; }
        std::size_t size() const { return end() - begin() ; }

        const T& operator[] ( std::size_t pos ) const { return begin()[pos] ; }
        const T& at( std::size_t pos ) const
        { if( pos >= size() ) throw std::out_of_range("matrix<T>::const_row: out of range element access") ; return begin()[pos] ; }

        const const_iterator begin_ ;
        const const_iterator end_ ;
    };

    row operator[] ( std::size_t pos ) { return { data_.begin() + pos*ncols, data_.begin() + (pos+1)*ncols } ; }
    const_row operator[] ( std::size_t pos ) const { return { data_.cbegin() + pos*ncols, data_.cbegin() + (pos+1)*ncols } ; }
    
    row at( std::size_t pos ) 
    {
        if( pos >= rows() ) throw std::out_of_range("matrix<T>::const_row: out of range element access") ;
        return { data_.begin() + pos*ncols, data_.begin() + (pos+1)*ncols } ; 
    }
    const_row at( std::size_t pos ) const 
    {
        if( pos >= rows() ) throw std::out_of_range("matrix<T>::const_row: out of range element access") ;
        return { data_.cbegin() + pos*ncols, data_.cbegin() + (pos+1)*ncols } ; 
    }

    // uBLAS style element access through call operator ()
    // http://www.boost.org/doc/libs/1_57_0/libs/numeric/ublas/doc/matrix.html#matrix
    // checked access; repace at() with operator[]() to change to unchecked access
    T& operator() ( std::size_t r, std::size_t c ) { return this->at(r).at(c) ; }
    const T& operator() ( std::size_t r, std::size_t c ) const { return this->at(r).at(c) ; }

    matrix() = default ;
    
    matrix( std::size_t nrows, std::size_t ncols, const T& v = T() ) : nrows(nrows), ncols(ncols), data_( nrows*ncols, v ) {}
    
    template < typename U > matrix( std::initializer_list< std::initializer_list<U> > ilist )
    {
        for( const auto& row : ilist )  ncols = std::max( ncols, row.size() ) ;
        data_.reserve( ilist.size() * ncols ) ;

        for( const auto& row : ilist )
        {
            ++nrows ;
            for( const U& v : row ) data_.emplace_back(v) ;
            data_.resize( nrows * ncols ) ; // fill up empty slots in the row
        }
    }

    matrix( const matrix<T>& that ) = default ; // defaulted copy constructor
    matrix( matrix<T>&& that ) noexcept = default ; // defaulted move constructor

    template < typename U > matrix( const matrix<U>& that,
                                    const typename std::enable_if< std::is_convertible<U,T>::value >::type* = nullptr )
                              : nrows( that.rows() ), ncols( that.cols() )
    {
        data_.resize( that.data().size() ) ;
        std::copy( std::begin( that.data() ), std::end( that.data() ), std::begin(data_) ) ;
    }

    matrix<T>& operator= ( const matrix<T>& that ) = default ; // defaulted copy assignment
    matrix<T>& operator= ( matrix<T>&& that ) noexcept = default ; // defaulted move assignment

    template < typename U >
    typename std::enable_if< std::is_convertible<U,T>::value, matrix<T>& >::type  operator= ( const matrix<U>& that )
    { return *this = matrix<T>(that) ; } // not the most efficient

    // implicitly declared destructor

    private:
        std::size_t nrows = 0 ;
        std::size_t ncols = 0 ;
        std::vector<T> data_ ;
};

template < typename C, typename TRAITS, typename T >
std::basic_ostream<C,TRAITS>& operator<< ( std::basic_ostream<C,TRAITS>& stm, const matrix<T>& m )
{
   for( std::size_t i = 0 ; i < m.rows() ; ++i )
   {
       for( const T& v : m[i] ) stm << v << stm.widen(' ') ;
       stm << stm.widen('\n') ;
   }
   return stm ;
}

/////////////////////////////////////////////////
// explicit instantiations for testing purposes
template struct matrix<long long> ;
template struct matrix<float>::row ;
template struct matrix<int>::const_row ;
//////////////////////////////////////////////////

int main()
{
    matrix<int> m1 { { {0,1,2,3}, {4,5,6,7} } };

    matrix<int> m2(m1) ;
    std::cout << m2 << '\n' ;

    matrix<double> m3(m1) ;
    std::cout << std::fixed << std::setprecision(1) << m3 << '\n' ;

    m2 = matrix<int> { { {0,1,2}, {3}, {6,7,8} } };
    std::cout << m2 << '\n' ;

    m3 = m2 ;
    std::cout << m3 << '\n' ;

    m2[1][1] = m2[1][2] = 9 ;
    std::cout << m2 << '\n' ;
}

http://coliru.stacked-crooked.com/a/5d5ab731b87121cc
cmajor28 (48)
I finished my matrix class using c style operations and I was going to post it but sadly I can't post attachments. I decided to create a variation of the class using vectors and what's above is a great start! Thanks a lot!

And I think once I'm done, it should support templates of class types. For example, if someone has a fraction class, with std::vector it wouldn't cause a problem. My version only supports int, float, etc. and complex numbers. It might support others but problems will likely arise.

One quick question though.

Why is this not valid
1
2
3
4
5
6
7
8
9
10
11
12
13
14
 
template <class T>
Matrix<T> Matrix<T>::trans() const {
	
	// Allocate memory for array
	T* array = (T*)malloc(m*n*sizeof(T));
	
	// Set array elements to opposite index of matrix element
	for (size_t i=0;i<m;i++) {
		for (size_t j=0;j<n;j++) {
			*(array+n*i+j) = *(matrix+n*j+i);
		}
	}
	return Matrix<T>(array,m,n);
}


While this is
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
 
template <class T>
Matrix<T> Matrix<T>::trans() const {
	
	// Allocate memory for array
	T* array = (T*)malloc(m*n*sizeof(T));
	
	// Set array elements to opposite index of matrix element
	for (size_t i=0;i<m;i++) {
		for (size_t j=0;j<n;j++) {
			*(array+n*i+j) = *(matrix+n*j+i);
		}
	}
        Matrix<T> temp(array,m,n);
	return temp;
}
Last edited on
JLBorges (13770)
> Why is this not valid return Matrix<T>(array,m,n);
> While this is Matrix<T> temp(array,m,n); return temp;

Either the copy constructor of Matrix<T> is not const-correct
or the move constructor of Matrix<T> has been explicitly deleted
(or both).
prestokeys (357)
matrix( matrix<T>&& that ) /*noexcept*/ = default ; // defaulted move constructor

What's is GCC 4.8.1's complaint with noexcept?

[Error] function 'matrix<T>& matrix<T>::operator=(matrix<T>&&) [with T = int]' defaulted on its first declaration with an exception-specification that differs from the implicit declaration 'matrix<int>& matrix<int>::operator=(matrix<int>&&)'
JLBorges (13770)
Fixed in g++ 4.9 http://coliru.stacked-crooked.com/a/efe3748d8f74d776

Even this does not work with g++ 4.8 ( doing too much in phase one of the two-phase look up? ).
1
2
3
 
matrix( matrix<T>&& that ) noexcept( std::is_nothrow_move_constructible< std::vector<T> >::value ) = default ;

matrix<T>& operator= ( matrix<T>&& that ) noexcept( std::is_nothrow_move_assignable< std::vector<T> >::value ) = default ;

http://coliru.stacked-crooked.com/a/74714105fe7a13bd
Last edited on
cmajor28 (48)
Thanks JLBorges. That was it.
Topic archived. No new replies allowed.