C++テンプレート実戦8:マトリクス乗算


行列乗算は反復器で実現され、行列は行優先方式で格納され、その重要な操作は「行」である.×行は1つの反復器によって移動を完了し、列には1つの列反復器によって移動を完了し、乗算はtransformによって完了します.乗算は累積反復器によって完了する必要があります.
1行列乗算で行に関連×列、行列が行優先方式で格納される場合、行の移動は比較的簡単で、列の移動は比較的複雑で、列の移動に対して反復器を採用して実現し、以下のようにする.
//    :skip_iterator.hpp
#pragma once
#include 

template
class skip_iterator :
    public std::iterator<:random_access_iterator_tag t="">
{
    T *pos;                     //       ,pos             
    size_t step;                //   ,      ,                   
public:
    //     
    typedef std::iterator<:random_access_iterator_tag t=""> base_type;
    typedef typename base_type::difference_type difference_type;
    typedef typename base_type::reference reference;

    //     
    skip_iterator(T *pos, size_t step) : pos(pos), step(step) {}
    skip_iterator(const skip_iterator &i) : pos(i.pos), step(step) {}

    difference_type
    operator - (skip_iterator r) {return (pos - r.pos) / step;}

    skip_iterator
    operator + (typename base_type::difference_type n)
    {return skip_iterator(pos + n * step);}

    skip_iterator
    operator ++() {pos += step; return *this;}//                                    step

    bool operator == (skip_iterator const &r) {return pos == r.pos;}
    bool operator != (skip_iterator const &r) {return !(*this == r);}

    //    
    reference operator *() {return *pos;}
};
   
2マトリクス類テンプレート、マトリクスは行優先方式で記憶され、特殊な設計が必要なのは列の移動であり、以下の通りである.
//    :matrix.hpp
#pragma once
#include "skip_iterator.hpp"
#include 

template
class matrix
{
public:
    //       
    typedef T value_type;
    typedef T* iterator; //       
    typedef T* row_iterator; //       
    typedef skip_iterator col_iterator; //       ,           skip_iterator
    //          
    typedef const T* const_iterator;
    typedef const T* const_row_iterator;
    typedef skip_iterator const_col_iterator;

private:
    T *data; //     
    size_t n_row;   //   
    size_t n_col;   //   
public:
    //      
    matrix(size_t n_row, size_t n_col) :
        data(new T[n_row * n_col]), n_row(n_row), n_col(n_col) {}

    //       
    matrix(matrix const &m) :
        data(new T[m.n_row * m.n_col]),
        n_row(m.n_row), n_col(m.n_col) {
        std::copy(m.begin(), m.end(), begin());
    }

    template
    matrix(size_t n_row, size_t n_col, Iterator i) :
        data(new T[n_row * n_col]), n_row(n_row), n_col(n_col) 
    {
        Iterator j = i;
        std::advance(j, n_row * n_col);
        std::copy(i, j, begin());
        //       C++11  copy_n
    }
    ~matrix() {delete[] data;} //       data   

    
    iterator begin() {return data;}
    iterator end() {return data + n_row * n_col;}
    row_iterator row_begin(size_t n) {return data + n * n_col;}// n    
    row_iterator row_end(size_t n) {return row_begin(n) + n_col;}// n    
    col_iterator col_begin(size_t n) {return col_iterator(data + n, n_col);}// n    ,n_col step
    col_iterator col_end(size_t n) {return col_begin(n) + n_row;}

    const_iterator begin() const {return data;}//    
    const_iterator end() const {return data + n_row * n_col;}
    const_row_iterator row_begin(size_t n) const {return data + n * n_col;}
    const_row_iterator row_end(size_t n) const {return row_begin(n) + n_col;}
    const_col_iterator col_begin(size_t n) const {
        return const_col_iterator(data + n, n_col);
    }
    const_col_iterator col_end(size_t n) const {return col_begin(n) + n_row;}

    size_t num_row() const {return n_row;}
    size_t num_col() const {return n_col;}

    T& operator() (size_t i, size_t j) {return data[i * n_col + j];}
    T const& operator() (size_t i, size_t j) const {return data[i * n_col + j];}

    //             
    matrix&
    operator=(matrix const &m) {
        if (&m == this) return *this;
        if (n_row * n_col < m.n_row * m.n_col) {
            delete[] data;
            data = new T[m.n_row * m.n_col];
        }
        n_row = m.n_row;
        n_col = m.n_col;
        std::copy(m.begin(), m.end(), begin());
    }
};

3累積反復器:行列の行×列は、2つのシーケンスの要素ごとに乗算された再累積に相当し、次のように、累積操作を実行するには反復器を使用します.
//    :accumulate_iterator.hpp
#pragma once
#include 
template
class accumulate_iterator :
    public std::iterator<:output_iterator_tag t="">
{
    T &ref_x;                   //         
    BinFunc bin_func;           //       。ref_x = bin_func(ref_x, v)
public:
    accumulate_iterator(T &ref_x, BinFunc bin_func) :
        ref_x(ref_x), bin_func(bin_func) {}

    //          
    accumulate_iterator operator*() {return *this;}

    //         
    template
    T0 const & operator=(T0 const &v) {ref_x = bin_func(ref_x, v);}

    accumulate_iterator& operator++() {return *this;}
};

//   accumulate_iterator     
template
accumulate_iterator
accumulater(T &ref_x, BinFunc bin_func)
{
    return accumulate_iterator(ref_x, bin_func);
}

4行列乗算の実現は以下の通りである.
#include "matrix.hpp"
#include "accumulate_iterator.hpp"
#include 
#include 

template
matrix operator * (matrix const &m0, matrix const &m1)
    throw (std::runtime_error)
{
    //         ,    。         
    if (m0.num_col() != m1.num_row())
        throw std::runtime_error("Bad matrix size for multiplication.");

    matrix m(m0.num_row(), m1.num_col());

    typename matrix::iterator pos = m.begin();
    
    for (size_t i = 0; i < m.num_row(); ++i) {
        for (size_t j = 0; j < m.num_col(); ++j) {
            *pos = 0;
            std::transform(m0.row_begin(i), m0.row_end(i),m1.col_begin(j),accumulater(*pos, std::plus()),std::multiplies());//    transform   #1#
            ++pos;
        }
    }

    return m;
}
#1#で説明:
template 
  OutputIterator transform (InputIterator first1, InputIterator last1,InputIterator first2,OutputIterator result, BinaryOperator binary_op)
{
  while (first1 != last1) {
    *result=binary_op(*first1,*first2++);//    #1# result      ,          。first2      
    ++result;//             
    ++first1;
  }
  return result;
}