trying to build a basic linear algebra library

[E093605o]

I am trying to implement a function to generate a Matrix as it is defined in my Matrix.h file initialized with random values drawn from a gaussian distribution.
However, when I try to print a randomly generated Matrix like that I get only (-)inf values. Can somebody tell me why and how I could fix this? Maybe a completely different approach rather than the struct approach I take.

//Matrix.h
#pragma once
#include <vector>
#include <cmath>
#include <cassert>
#include <iostream>
#include <tuple>
#include <random>

template<typename Type>
class Matrix {

  size_t cols{};
  size_t rows{};

 public:
  std::vector<std::vector<Type>> data;
  //std::tuple<size_t, size_t> shape;
  size_t elementCount{};

  /* constructors */
  Matrix(size_t rowsArg, size_t colsArg) : cols(colsArg), rows(rowsArg)
  {
		data = std::vector<std::vector<Type>>(rows,std::vector<Type>(cols));
	}

  Matrix(){};

  //methods
  void print();

  Matrix<Type> matmul(Matrix<Type> &m);


  Type& operator()(size_t row, size_t col) {
    assert(row < data.size() && col < data[0].size());
    return data[row][col];
}
  

};

// methods
template<typename Type>
void Matrix<Type>::print(){
  for (int i = 0; i < rows; i++){
    for (int j = 0; j < cols; j++){
      std::cout << data[i][j] << " ";
    }
    std::cout << std::endl;
  }
}

template <typename Type>
Matrix<Type> Matrix<Type>::matmul(Matrix<Type> &target) {
    assert(cols == target.rows);
    Matrix output(rows, target.cols);

    for (size_t r = 0; r < output.rows; ++r) {
      for (size_t c = 0; c < output.cols; ++c) {
        for (size_t k = 0; k < target.rows; ++k)
          output(r, c) += (*this)(r, k) * target(k, c);
      }
    }
    return output;
  };

template <typename T>
struct mtx {
  static Matrix<T> randn(size_t rows, size_t cols) {
    Matrix<T> M(rows, cols);

    std::random_device rd{};
    std::mt19937 gen{rd()};

    // init Gaussian distr. w/ N(mean=0, stdev=1/sqrt(numel))
    T n(M.elementCount);
    T stdev{1 / sqrt(n)};
    std::normal_distribution<T> d{0, stdev};

    // fill each element w/ draw from distribution
    for (size_t r = 0; r < rows; ++r) {
      for (int c = 0; c < cols; ++c) {
        M(r, c) = d(gen);
      }
    }
    return M;
  }
};

//Matrix.cpp
#include "Matrix.h"

int main(){
  auto A = mtx<double>::randn(2,3);
  auto B = mtx<double>::randn(3,2);

  // all of the print() calls print (-)inf values or NaN for A.matmul(B).print();
  A.print();
  B.print();
  A.matmul(B).print();
  return 0;
}

[seeplus]

You don't set elementCount!

As 1 file, then perhaps:

#include <vector>
#include <cmath>
#include <cassert>
#include <iostream>
#include <random>

template<typename Type>
class Matrix {
	std::vector<std::vector<Type>> data;

public:
	const size_t cols {};
	const size_t rows {};
	const size_t elementCount {};

	Matrix(size_t rowsArg, size_t colsArg) : cols(colsArg), rows(rowsArg), elementCount(rowsArg * colsArg),
		data(std::vector<std::vector<Type>>(rowsArg, std::vector<Type>(colsArg))) {}

	Matrix() {}

	void print() const;

	Matrix<Type> matmul(const Matrix<Type>& m) const;

	Type& operator()(size_t row, size_t col) {
		assert(row < data.size() && col < data[0].size());
		return data[row][col];
	}

	Type operator()(size_t row, size_t col) const {
		assert(row < data.size() && col < data[0].size());
		return data[row][col];
	}
};

template<typename Type>
void Matrix<Type>::print() const {
	for (const auto& r : data) {
		for (const auto& e : r)
			std::cout << e << ' ';

		std::cout << '\n';
	}

	std::cout << '\n';
}

template <typename Type>
Matrix<Type> Matrix<Type>::matmul(const Matrix<Type>& target) const {
	assert(cols == target.rows);
	Matrix output(rows, target.cols);

	for (size_t r {}; r < output.rows; ++r) {
		for (size_t c {}; c < output.cols; ++c)
			for (size_t k {}; k < target.rows; ++k)
				output(r, c) += (r, k) * target(k, c);
	}

	return output;
};

template <typename T>
struct mtx {
	static inline std::mt19937 gen { std::random_device {}() };

	static Matrix<T> randn(size_t rows, size_t cols) {
		Matrix<T> M(rows, cols);
		const T stdev { 1.0 / std::sqrt(static_cast<T>(M.elementCount)) };
		std::normal_distribution<T> d { 0, stdev };

		for (size_t r {}; r < rows; ++r)
			for (size_t c {}; c < cols; ++c)
				M(r, c) = d(gen);

		return M;
	}
};

int main() {
	const auto A { mtx<double>::randn(2, 3) };
	const auto B { mtx<double>::randn(3, 2) };

	A.print();
	B.print();

	A.matmul(B).print();
}

0.0285513 0.193599 -0.919214
0.321249 0.0546027 -0.154442

-0.187649 0.431069
1.21961 0.246825
-0.227973 -0.679934

0.440315 0.685098
0.0415206 0.256968

[JLBorges]

Why do we need to store these? (They can be inexpensively computed.)

template<typename Type>
class Matrix {
     // ...
	const size_t cols {};
	const size_t rows {};
	const size_t elementCount {};
     // ...
};

For example:

#include <vector>
#include <cmath>
#include <cassert>
#include <iostream>
#include <iomanip>
#include <cassert>
#include <random>
#include <cstdlib>

template < typename T > class Matrix
{
    public:

        using row_type = std::vector<T> ;

        Matrix( std::size_t nrows, std::size_t ncols, const T& v = {} ) : data( nrows, row_type( ncols, v ) ) {}

        std::size_t nrows() const noexcept { return data.size() ; }
        std::size_t ncols() const noexcept { return data.empty() ? 0 : data.front().size() ; }
        std::size_t element_count() const noexcept { return nrows() * ncols() ; }

        const row_type& operator[] ( std::size_t p ) const { return data[p] ; }
        const row_type& at( std::size_t p ) const { return data.at(p) ; }

        template < typename UNARY_FN > void for_each( UNARY_FN fn ) const
        { for( const row_type& row : data ) for( const T& v : row ) fn(v) ; }

        template < typename NULLARY_FN > void generate( NULLARY_FN fn )
        { for( row_type& row : data ) for( T& v : row ) v = fn() ; }

        void fill( const T& v )
        { for( row_type& row : data ) for( T& elem : row ) elem = v ; }

        template < typename UNARY_FN > void transform( UNARY_FN fn )
        { for( row_type& row : data ) for( T& v : row ) v = fn(v) ; }

        // etc.

        bool valid() const noexcept // class invariant. TO DO: add assertions
        {
            for( const row_type& row : data ) if( row.size() != ncols() ) return false ;
            return true ;
        }

    private: std::vector<row_type> data ;

    friend std::ostream& operator<< ( std::ostream& stm, const Matrix<T>& mtx )
    {
        for( const auto& row : mtx.data )
        {
            // TO DO: compute appropriate width at run-time
            for( const T& v : row ) stm << std::setw(8) << v ;
            stm << '\n' ;
        }
        return stm ;
    }

    // etc.
};

int main()
{
    Matrix<int> mtx( 6, 9, 22 ) ;
    assert( mtx.valid() && mtx.nrows() == 6 && mtx.ncols() == 9 && mtx.element_count() == 6*9 ) ;
    for( std::size_t r = 0 ; r < mtx.nrows() ; ++r )
        for( std::size_t c = 0 ; c < mtx.ncols() ; ++c )
            assert( mtx[r][c] == 22 ) ;
    std::cout << mtx << "\n---------------------\n" ;

    std::mt19937 rng( std::random_device{}() ) ;
    mtx.generate( [&] { return rng()%100 - 50 ; } ) ;
    std::cout << mtx << "\n---------------------\n" ;

    mtx.transform( []( int v ) { return std::abs(v) % 10 ; } ) ;
    std::cout << mtx << "\n---------------------\n" ;
}

https://coliru.stacked-crooked.com/a/11c9017c09138770

[E093605o]

@seeplus thanks I didnt see that as usual...

@JLBorges maybe - I am 3 days into learning c++ so Ill try to focus on getting a program to behave as intended for now.