// Diggins Public Domain Post #5
// Contributed to Public Domain by Christopher Diggins
// This is an efficient matrix class implementation for when rows and
columns are known at compile-time
#include <algorithm>
#include <iterator>
#include <cassert>
//#include "static_assert.hpp"
template<class Value_T, int Size_N, int Step_N>
class kslice_iter
{
public:
// public typedefs
typedef Value_T value_type;
typedef Value_T& reference;
typedef const Value_T& const_reference;
typedef std::size_t size_type;
typedef std::ptrdiff_t difference_type;
typedef kslice_iter self;
// constructors
explicit kslice_iter(Value_T* x) : m(x) {
}
// operators
self& operator++() {
m += step();
return *this;
}
self operator++(int) {
self tmp = *this;
m += step();
return tmp;
}
Value_T& operator[](size_type n) {
assert(n < size());
return m[n * Step_N];
}
const Value_T& operator[](size_type n) const {
assert(n < N);
return m[n * Step_N];
}
static int size() {
return Size_N;
}
static int step() {
return Step_N;
}
Value_T& operator*() {
return *m;
}
friend difference_type operator-(const self& x, const self& y) {
return (y.m - x.m) / step();
}
friend bool operator==(const self& x, const self& y) {
return x.m == y.m;
}
friend bool operator!=(const self& x, const self& y) {
return x.m != y.m;
}
friend bool operator<(const self& x, const self& y) {
return x.m < y.m;
}
private:
// default constructor hidden from view
kslice_iter() { };
Value_T* m;
};
template<int N, class T, class RowIter_T, class ColIter_T>
T dot_product(RowIter_T x, ColIter_T y)
{
T ret(0);
for (int i=N; i > 0; --i) {
ret += *x++ * *y++;
}
return ret;
}
template<class T, class M1, class M2, class M3>
kmatrix_multiply(const M1& m1, const M2& m2, M3& result)
{
//STATIC_ASSERT(M1::cols == M2::rows);
//STATIC_ASSERT(M3::rows == M1::rows);
//STATIC_ASSERT(M3::cols == M2::cols);
for (int i=0; i < M1::rows; ++i) {
for (int j=0; j < M2::cols; ++j) {
result[i][j] = dot_product<M1::cols, T>(m1.row(i), m2.column(j));
}
}
}
template<class Value_T, int Rows_N, int Cols_N>
class kmatrix
{
public:
// public typedefs
typedef Value_T value_type;
typedef kmatrix self;
typedef Value_T* iterator;
typedef const Value_T* const_iterator;
typedef kslice_iter<Value_T, Cols_N, 1> row_type;
typedef kslice_iter<Value_T, Rows_N, Cols_N> col_type;
// public constants
static const int rows = Rows_N;
static const int cols = Cols_N;
// constructors
kmatrix() {
std::fill(begin(), end(), Value_T(0));
}
kmatrix(Value_T& x) {
std::fill(begin(), end(), x);
}
kmatrix(const self& x) {
operator=(x);
}
// public functions
self& operator=(const self& x) {
std::copy(x.begin(), x.end(), begin());
return *this;
}
self& operator=(Value_T& x) {
std::fill(begin(), end(), x);
return *this;
}
row_type row(int n) const {
assert(n < rows);
return row_type(m + (n * Cols_N));
}
col_type column(int n) const {
assert(n < cols);
return col_type(m + n);
}
row_type operator[](int n) {
return row(n);
}
const row_type operator[](int n) const {
return row(n);
}
iterator begin() {
return m;
}
iterator end() {
return m + size();
}
const_iterator begin() const {
return m;
}
const_iterator end() const {
return m + size();
}
static int size() {
return Rows_N * Cols_N;
}
// operators
self& operator+=(const self& x) {
for (int i=0; i < size(); ++i) {
m[i] += x.m[i];
}
return *this;
}
self& operator-=(const self& x) {
for (int i=0; i < size(); ++i) {
m[i] -= x.m[i];
}
return *this;
}
self& operator+=(value_type x) {
for (int i=0; i < size(); ++i) {
m[i] += x;
}
return *this;
}
self& operator-=(value_type x) {
for (int i=0; i < size(); ++i) {
m[i] -= x;
}
return *this;
}
self& operator*=(value_type x) {
for (int i=0; i < size(); ++i) {
m[i] *= x;
}
return *this;
}
self& operator/=(value_type x) {
for (int i=0; i < size(); ++i) {
m[i] /= x;
}
return *this;
}
self& operator=(value_type x) {
std::fill(begin(), end(), x);
return *this;
}
self operator-() {
self x;
for (int i=0; i < size(); ++i) {
x.m[i] = -m[i];
}
return x;
}
// friends
friend self operator+(const self& x, const self& y) {
return self(x) += y;
}
friend self operator-(const self& x, const self& y) {
return self(x) -= y;
}
friend self operator+(const self& x, value_type y) {
return self(x) += y;
}
friend self operator-(const self& x, value_type y) {
return self(x) -= y;
}
friend self operator*(const self& x, value_type y) {
return self(x) *= y;
}
friend self operator/(const self& x, value_type y) {
return self(x) /= y;
}
template<int Cols2_N>
friend kmatrix<Value_T, Rows_N, Cols2_N>
operator*(const self& x, const kmatrix<Value_T, Cols_N, Cols2_N>& y) {
kmatrix<Value_T, Rows_N, Cols2_N> ret;
kmatrix_multiply(x, y, ret);
return ret;
}
private:
mutable Value_T m[Rows_N * Cols_N];
}; 1 1370
// public domain matrix benchmark code
#include "kmatrix.hpp"
#include <cassert>
#include <iostream>
#include <cstdlib>
#include <boost/numeric/ublas/matrix.hpp>
class timer
{
public:
timer() : start(std::clock()) {
}
double elapsed() const {
return (double(std::clock() - start) / CLOCKS_PER_SEC) * 1000;
}
private:
std::clock_t start;
};
void output_time_elapsed(const timer& t, std::ostream& o = std::cout) {
o << static_cast<int>(t.elapsed()) << " msec elapsed" << std::endl;
}
struct reporting_timer {
~reporting_timer() {
output_time_elapsed(m);
}
timer m;
};
template<int Rows_N, int Cols_N, class Matrix_T>
void init_kmatrix(Matrix_T& m)
{
for (int i=0; i < Rows_N; ++i) {
for (int j=0; j < Cols_N; ++j) {
m[i][j] = i * Cols_N + j;
}
}
}
template<int Rows_N, int Cols_N, class Matrix_T>
void init_ublas_matrix(Matrix_T& m)
{
for (int i=0; i < Rows_N; ++i) {
for (int j=0; j < Cols_N; ++j) {
m(i, j) = i * Cols_N + j;
}
}
}
template<class T, int Iters_N, int M, int R, int N>
void kmatrix_performance_test()
{
std::cout << "kmatrix ";
kmatrix<T, M, R> km1;
kmatrix<T, R, N> km2;
kmatrix<T, M, N> km3;
init_kmatrix<M, R>(km1);
init_kmatrix<R, M>(km2);
{
reporting_timer t;
for (int i=0; i < Iters_N; ++i)
{
kmatrix_multiply<T>(km1, km2, km3);
}
}
std::cout << "ublas ";
boost::numeric::ublas::matrix<T> m1(M,R);
boost::numeric::ublas::matrix<T> m2(R,N);
boost::numeric::ublas::matrix<T> m3(M,N);
init_ublas_matrix<M, R>(m1);
init_ublas_matrix<R, M>(m2);
{
reporting_timer t;
for (int i=0; i < Iters_N; ++i)
{
boost::numeric::ublas::noalias(m3) =
boost::numeric::ublas::prod(m1,m2);
}
}
for (int i=0; i<M; ++i) {
for (int j=0; j<N; ++j) {
T x = m3(i,j);
T y = km3[i][j];
ASSURE(x == y);
}
}
}
void benchmarks()
{
std::cout << "comparing integer matricies" << std::endl;
std::cout << "2,2 X 2,2" << std::endl;
kmatrix_performance_test<int, 100000, 2, 2, 2>();
std::cout << "3,3 X 3,3" << std::endl;
kmatrix_performance_test<int, 100000, 3, 3, 3>();
std::cout << "100,100 X 100,100" << std::endl;
kmatrix_performance_test<int, 10, 100, 100, 100>();
std::cout << "100,1 X 1,100" << std::endl;
kmatrix_performance_test<int, 100, 100, 1, 100>();
std::cout << "1,100 X 100, 1" << std::endl;
kmatrix_performance_test<int, 100000, 1, 100, 1>();
std::cout << "comparing double matricies" << std::endl;
std::cout << "2,2 X 2,2" << std::endl;
kmatrix_performance_test<double, 100000, 2, 2, 2>();
std::cout << "3,3 X 3,3" << std::endl;
kmatrix_performance_test<double, 100000, 3, 3, 3>();
std::cout << "100,100 X 100,100" << std::endl;
kmatrix_performance_test<double, 10, 100, 100, 100>();
std::cout << "100,1 X 1,100" << std::endl;
kmatrix_performance_test<double, 100, 100, 1, 100>();
std::cout << "1,100 X 100, 1" << std::endl;
kmatrix_performance_test<double, 100000, 1, 100, 1>();
}
int main()
{
benchmarks();
system("pause");
return 0;
}
/*
my results:
comparing integer matricies
2,2 X 2,2
kmatrix 78 msec elapsed
ublas 172 msec elapsed
3,3 X 3,3
kmatrix 219 msec elapsed
ublas 437 msec elapsed
100,100 X 100,100
kmatrix 516 msec elapsed
ublas 1156 msec elapsed
100,1 X 1,100
kmatrix 125 msec elapsed
ublas 203 msec elapsed
1,100 X 100, 1
kmatrix 500 msec elapsed
ublas 1172 msec elapsed
comparing double matricies
2,2 X 2,2
kmatrix 79 msec elapsed
ublas 171 msec elapsed
3,3 X 3,3
kmatrix 250 msec elapsed
ublas 422 msec elapsed
100,100 X 100,100
kmatrix 656 msec elapsed
ublas 1344 msec elapsed
100,1 X 1,100
kmatrix 156 msec elapsed
ublas 219 msec elapsed
1,100 X 100, 1
kmatrix 703 msec elapsed
ublas 1172 msec elapsed
Press any key to continue . . .
*/ This thread has been closed and replies have been disabled. Please start a new discussion. Similar topics
by: SoloCDM |
last post by:
How do I keep my entire web page at a fixed width?
*********************************************************************
Signed,
SoloCDM
|
by: Mason A. Clark |
last post by:
Masters:
On two or three-column layouts, one column often has a list
of links. Scrolling the page hides them. I'm aware there's
supposed to be the ability to fix the column (frame-like).
I...
|
by: christopher diggins |
last post by:
Welcome to the first installment of the Diggins Public Domain Posts (PDP).
There is a significant dearth of good public domain C++ code. All of it
seems to come with one kind of a license or...
|
by: christopher diggins |
last post by:
// meta_binder.hpp
// The Diggins PDP (Public Domain Post) #2
// Public Domain code by Christopher Diggins, May 22, 2005
//
// Description:
// A meta-binder binds function objects to function...
|
by: christopher diggins |
last post by:
// binary_aritmetic.hpp
// Diggins PDP (Public Domain Post) #1.2
// by Christopher Diggins, May 24, 2005
//
// Comment:
// Provides several algorithms for unsigned binary arithmetic
//
//...
|
by: christopher diggins |
last post by:
// big_int.hpp
// The Diggins PDP (Public Domain Post) #3
// Public Domain code by Christopher Diggins, May 22, 2005
//
// Description:
// A naively implemented unsigned integer type for...
|
by: christopher diggins |
last post by:
// The Diggins PDP (Public Domain Post) #4
// public domain code for computing the FFT
// contributed by Christopher Diggins, 2005
template<int N>
unsigned int bit_reverse(unsigned int x) {...
|
by: Otie |
last post by:
Hello,
I am using the MSFlexGrd Control in VB5. I have 1 fixed row and one fixed
column. I am trying to do a sort when the user clicks a column in the FIXED
ROW. But when I capture the row...
|
by: Jeff |
last post by:
Hey
I'm wondering how the Fixed-Width Text Format is
What I know is that the top line in this text format will contain column
names. and each row beneath the top line represent for example a...
|
by: isladogs |
last post by:
The next Access Europe meeting will be on Wednesday 6 Mar 2024 starting at 18:00 UK time (6PM UTC) and finishing at about 19:15 (7.15PM).
In this month's session, we are pleased to welcome back...
|
by: Vimpel783 |
last post by:
Hello!
Guys, I found this code on the Internet, but I need to modify it a little. It works well, the problem is this: Data is sent from only one cell, in this case B5, but it is necessary that data...
|
by: jfyes |
last post by:
As a hardware engineer, after seeing that CEIWEI recently released a new tool for Modbus RTU Over TCP/UDP filtering and monitoring, I actively went to its official website to take a look. It turned...
|
by: ArrayDB |
last post by:
The error message I've encountered is; ERROR:root:Error generating model response: exception: access violation writing 0x0000000000005140, which seems to be indicative of an access violation...
|
by: CloudSolutions |
last post by:
Introduction:
For many beginners and individual users, requiring a credit card and email registration may pose a barrier when starting to use cloud servers. However, some cloud server providers now...
|
by: Defcon1945 |
last post by:
I'm trying to learn Python using Pycharm but import shutil doesn't work
|
by: Shællîpôpï 09 |
last post by:
If u are using a keypad phone, how do u turn on JavaScript, to access features like WhatsApp, Facebook, Instagram....
|
by: Faith0G |
last post by:
I am starting a new it consulting business and it's been a while since I setup a new website. Is wordpress still the best web based software for hosting a 5 page website? The webpages will be...
|
by: isladogs |
last post by:
The next Access Europe User Group meeting will be on Wednesday 3 Apr 2024 starting at 18:00 UK time (6PM UTC+1) and finishing by 19:30 (7.30PM).
In this session, we are pleased to welcome former...
| |