Lösungsvorschlag

Zur ersten Frage: Das RandomValues-Objekt ist zunächst ein temporäres Objekt, das dann durch die Parameterübergabe kopiert und dadurch in initGeMatrix zu einem lokalen Objekt wird. Bei der Parameterübergabe wird implizit der copy constructor aufgerufen. Da kein expliziter gegeben ist, werden per Voreinstellungen alle einzelnen Komponenten, also mt und uniform, kopiert.

Zur zweiten Frage: Der Datentyp std::random_device hat keinen copy constructor. Deswegen klappt die Parameterübergabe nicht, wenn es in RandomValues es eine entsprechend deklarierte Variable gibt. Deswegen wurde hier std::random_device()() verwendet. Das erste Klammernpaar ist für den Aufruf des Konstruktors. Bei dem so temporär erzeugten Objekt wird dann der Funktionsoperator aufgerufen. Hierzu ist das zweite Klammernpaar notwendig.

Beim dritten Punkt war es nur notwendig, die Funktion initGeMatrix in bench.h zu vereinfachen, d.h. der dritte und vierte Parameter wurden einfach gestrichen:

template <typename T, typename Index, typename InitValue>
void
initGeMatrix(Index m, Index n, T *A, Index incRowA, Index incColA,
             InitValue initValue)
{
    for (Index j=0; j<n; ++j) {
        for (Index i=0; i<m; ++i) {
            A[i*incRowA+j*incColA] = initValue(i, j);
        }
    }
}

So lässt sich die RandomValues-Klasse in Bezug auf die vierte Frage um einen weiteren Konstruktor für reproduzierbare Zahlenfolgen ergänzen:

template<typename T, typename Index>
struct RandomValues {
   std::mt19937 mt;
   std::uniform_real_distribution<T> uniform;

   RandomValues() :
      mt(std::random_device()()), uniform(-100, 100) {
   }
   RandomValues(unsigned int seed) :
      mt(seed), uniform(-100, 100) {
   }
   T operator()(Index i, Index j) {
      return uniform(mt);
   }
};

Und so könnte IncreasingRandomValues aussehen:

template<typename T, typename Index>
struct IncreasingRandomValues {
   std::mt19937 mt;
   std::uniform_real_distribution<T> uniform;
   T last_value;

   IncreasingRandomValues() :
      mt(std::random_device()()), uniform(0, 1<<16), last_value(0) {
   }
   IncreasingRandomValues(unsigned int seed) :
      mt(seed), uniform(1, 1<<16) {
   }
   T operator()(Index i, Index j) {
      last_value += uniform(mt);
      return last_value;
   }
};

Hier sind alle Dateien zum Lösungsvorschlag:

session11/step04/bench.h

#ifndef HPC_BENCH_H
#define HPC_BENCH_H 1

#include <chrono>
#include <cmath>
#include <complex>
#include <random>

namespace bench {

template <typename T, typename Index>
T
asumDiffGeMatrix(Index m, Index n,
                 const T *A, Index incRowA, Index incColA,
                 T *B, Index incRowB, Index incColB)
{

    T asum = 0;

    for (Index j=0; j<n; ++j) {
        for (Index i=0; i<m; ++i) {
            asum += std::abs(B[i*incRowB+j*incColB] - A[i*incRowA+j*incColA]);
        }
    }
    return asum;
}

template <typename T, typename Index, typename InitValue>
void
initGeMatrix(Index m, Index n, T *A, Index incRowA, Index incColA,
             InitValue initValue)
{
    for (Index j=0; j<n; ++j) {
        for (Index i=0; i<m; ++i) {
            A[i*incRowA+j*incColA] = initValue(i, j);
        }
    }
}

template <typename T>
struct WallTime
{
    void
    tic()
    {
        t0 = std::chrono::high_resolution_clock::now();
    }

    T
    toc()
    {
        using namespace std::chrono;

        elapsed = high_resolution_clock::now() - t0;
        return duration<T,seconds::period>(elapsed).count();
    }

    std::chrono::high_resolution_clock::time_point t0;
    std::chrono::high_resolution_clock::duration   elapsed;
};


} // namespace bench


#endif // HPC_BENCH_H

session11/step04/gematrix.h

#ifndef HPC_GEMATRIX_H
#define HPC_GEMATRIX_H 1

#include <cstddef>
#include <cassert>
#include "bench.h"

namespace matvec {

enum class StorageOrder {
    ColMajor,
    RowMajor
};

template <typename T, typename I>
    struct GeMatrixView;

template <typename T, typename I=std::size_t>
struct GeMatrix
{
    typedef T                       ElementType;
    typedef I                       Index;
    typedef GeMatrix<T,Index>       NoView;
    typedef GeMatrixView<T,Index>   View;

    GeMatrix(Index m, Index n, StorageOrder order=StorageOrder::ColMajor)
        : m(m), n(n),
          incRow(order==StorageOrder::ColMajor ? 1: n),
          incCol(order==StorageOrder::RowMajor ? 1: m),
          data(new T[m*n])
    {
    }

    ~GeMatrix()
    {
        delete[] data;
    }

    const ElementType &
    operator()(Index i, Index j) const
    {
        assert(i<m && j<n);
        return data[i*incRow + j*incCol];
    }

    ElementType &
    operator()(Index i, Index j)
    {
        assert(i<m && j<n);
        return data[i*incRow + j*incCol];
    }

    View
    operator()(Index i, Index j, Index numRows, Index numCols)
    {
        assert(i+numRows<=m);
        assert(j+numCols<=n);
        return View(numRows, numCols, &(operator()(i,j)), incRow, incCol);
    }

    const Index     m, n, incRow, incCol;
    ElementType*    data;
};

template <typename T, typename I=std::size_t>
struct GeMatrixView
{
    typedef T                       ElementType;
    typedef I                       Index;
    typedef GeMatrix<T,Index>       NoView;
    typedef GeMatrixView<T,Index>   View;

    GeMatrixView(Index m, Index n, T *data, Index incRow, Index incCol)
        : m(m), n(n), incRow(incRow), incCol(incCol), data(data)
    {
    }

    GeMatrixView(const GeMatrixView &rhs)
        : m(rhs.m), n(rhs.n), incRow(rhs.incRow), incCol(rhs.incCol),
          data(rhs.data)
    {
    }

    const ElementType &
    operator()(Index i, Index j) const
    {
        assert(i<m && j<n);
        return data[i*incRow + j*incCol];
    }

    ElementType &
    operator()(Index i, Index j)
    {
        assert(i<m && j<n);
        return data[i*incRow + j*incCol];
    }

    View
    operator()(Index i, Index j, Index numRows, Index numCols)
    {
        assert(i+numRows<=m);
        assert(j+numCols<=n);
        return View(numRows, numCols, &(operator()(i,j)), incRow, incCol);
    }

    const Index     m, n, incRow, incCol;
    ElementType*    data;
};

//
//  Interface for bench
//
template <typename GeMatrix, typename InitValue>
void
initGeMatrix(GeMatrix &A, InitValue initValue)
{
    bench::initGeMatrix(A.m, A.n, A.data, A.incRow, A.incCol, initValue);
}

template <typename GeMatrix>
typename GeMatrix::ElementType
asumDiffGeMatrix(const GeMatrix &A, const GeMatrix &B)
{
    return bench::asumDiffGeMatrix(A.m, A.n,
                                   A.data, A.incRow, A.incCol,
                                   B.data, B.incRow, B.incCol);
}

} // namespace matvec

#endif // HPC_GEMATRIX_H

session11/step04/test_initmatrix.cpp

#include <cstdlib>
#include <cstdio>
#include "gematrix.h"

template<typename T, typename Index>
T init_value(Index i, Index j, Index m, Index n) {
   return T(j)*T(n) + T(i) + T(1);
}

template<typename T, typename Index>
std::complex<T> init_complex_value(Index i, Index j, Index m, Index n) {
   return std::complex<T>(T(i), T(j));
}

void print_value(double value) {
   std::printf(" %10.1lf", value);
}

void print_value(std::complex<double> value) {
   std::printf(" (%4.1lf, %4.1lf)", value.real(), value.imag());
}

template<typename Matrix>
void print_matrix(const Matrix& A) {
   for (std::size_t i = 0; i < A.m; ++i) {
      std::printf("  ");
      for (std::size_t j = 0; j < A.n; ++j) {
         print_value(A(i, j));
      }
      std::printf("\n");
   }
}

template<typename T, typename Index>
struct RandomValues {
   std::mt19937 mt;
   std::uniform_real_distribution<T> uniform;

   RandomValues() :
      mt(std::random_device()()), uniform(-100, 100) {
   }
   RandomValues(unsigned int seed) :
      mt(seed), uniform(-100, 100) {
   }
   T operator()(Index i, Index j) {
      return uniform(mt);
   }
};

template<typename T, typename Index>
struct IncreasingRandomValues {
   std::mt19937 mt;
   std::uniform_real_distribution<T> uniform;
   T last_value;

   IncreasingRandomValues() :
      mt(std::random_device()()), uniform(0, 1<<16), last_value(0) {
   }
   IncreasingRandomValues(unsigned int seed) :
      mt(seed), uniform(1, 1<<16) {
   }
   T operator()(Index i, Index j) {
      last_value += uniform(mt);
      return last_value;
   }
};

int main() {
   using namespace matvec;
   GeMatrix<double> A(3, 7, StorageOrder::ColMajor);
   initGeMatrix(A, RandomValues<double, std::size_t>(42));
   printf("A:\n");
   print_matrix(A);

   GeMatrix<double> B(2, 8, StorageOrder::ColMajor);
   initGeMatrix(B, RandomValues<double, std::size_t>(42));
   printf("B:\n");
   print_matrix(B);

   GeMatrix<double> C(3, 7, StorageOrder::ColMajor);
   initGeMatrix(C, IncreasingRandomValues<double, std::size_t>());
   printf("C:\n");
   print_matrix(C);
}