diff options
author | Gael Guennebaud <g.gael@free.fr> | 2009-02-10 18:06:05 +0000 |
---|---|---|
committer | Gael Guennebaud <g.gael@free.fr> | 2009-02-10 18:06:05 +0000 |
commit | cbbc6d940bd6ddf7d96352b4b4e3033ff4e555fe (patch) | |
tree | c9e905f59b380bac934bf44dd635060cb1f2c264 /test | |
parent | a0cc5fba0a3dc600e858720c9fcde3a9c13e740a (diff) |
* add ei_predux_mul internal function
* apply Ricard Marxer's prod() patch with fixes for the vectorized path
Diffstat (limited to 'test')
-rw-r--r-- | test/CMakeLists.txt | 1 | ||||
-rw-r--r-- | test/packetmath.cpp | 5 | ||||
-rw-r--r-- | test/prod.cpp | 86 |
3 files changed, 92 insertions, 0 deletions
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 4d2ac7501..a20a6a12f 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -149,6 +149,7 @@ ei_add_test(sparse_basic) ei_add_test(sparse_product) ei_add_test(sparse_solvers " " "${SPARSE_LIBS}") ei_add_test(reverse) +ei_add_test(prod) # print a summary of the different options message("************************************************************") diff --git a/test/packetmath.cpp b/test/packetmath.cpp index d746a350e..0397d9b28 100644 --- a/test/packetmath.cpp +++ b/test/packetmath.cpp @@ -124,6 +124,11 @@ template<typename Scalar> void packetmath() for (int i=0; i<PacketSize; ++i) ref[0] += data1[i]; VERIFY(ei_isApprox(ref[0], ei_predux(ei_pload(data1))) && "ei_predux"); + + ref[0] = 1; + for (int i=0; i<PacketSize; ++i) + ref[0] *= data1[i]; + VERIFY(ei_isApprox(ref[0], ei_predux_mul(ei_pload(data1))) && "ei_predux_mul"); for (int j=0; j<PacketSize; ++j) { diff --git a/test/prod.cpp b/test/prod.cpp new file mode 100644 index 000000000..12f8c0c48 --- /dev/null +++ b/test/prod.cpp @@ -0,0 +1,86 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2008 Benoit Jacob <jacob.benoit.1@gmail.com> +// +// Eigen is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 3 of the License, or (at your option) any later version. +// +// Alternatively, you can redistribute it and/or +// modify it under the terms of the GNU General Public License as +// published by the Free Software Foundation; either version 2 of +// the License, or (at your option) any later version. +// +// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS +// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License and a copy of the GNU General Public License along with +// Eigen. If not, see <http://www.gnu.org/licenses/>. + +#include "main.h" + +template<typename MatrixType> void matrixProd(const MatrixType& m) +{ + typedef typename MatrixType::Scalar Scalar; + + int rows = m.rows(); + int cols = m.cols(); + + MatrixType m1 = MatrixType::Random(rows, cols); + + VERIFY_IS_MUCH_SMALLER_THAN(MatrixType::Zero(rows, cols).prod(), Scalar(1)); + VERIFY_IS_APPROX(MatrixType::Ones(rows, cols).prod(), Scalar(1)); + Scalar x = Scalar(1); + for(int i = 0; i < rows; i++) for(int j = 0; j < cols; j++) x *= m1(i,j); + VERIFY_IS_APPROX(m1.prod(), x); +} + +template<typename VectorType> void vectorProd(const VectorType& w) +{ + typedef typename VectorType::Scalar Scalar; + int size = w.size(); + + VectorType v = VectorType::Random(size); + for(int i = 1; i < size; i++) + { + Scalar s = Scalar(1); + for(int j = 0; j < i; j++) s *= v[j]; + VERIFY_IS_APPROX(s, v.start(i).prod()); + } + + for(int i = 0; i < size-1; i++) + { + Scalar s = Scalar(1); + for(int j = i; j < size; j++) s *= v[j]; + VERIFY_IS_APPROX(s, v.end(size-i).prod()); + } + + for(int i = 0; i < size/2; i++) + { + Scalar s = Scalar(1); + for(int j = i; j < size-i; j++) s *= v[j]; + VERIFY_IS_APPROX(s, v.segment(i, size-2*i).prod()); + } +} + +void test_prod() +{ + for(int i = 0; i < g_repeat; i++) { + CALL_SUBTEST( matrixProd(Matrix<float, 1, 1>()) ); + CALL_SUBTEST( matrixProd(Matrix2f()) ); + CALL_SUBTEST( matrixProd(Matrix4d()) ); + CALL_SUBTEST( matrixProd(MatrixXcf(3, 3)) ); + CALL_SUBTEST( matrixProd(MatrixXf(8, 12)) ); + CALL_SUBTEST( matrixProd(MatrixXi(8, 12)) ); + } + for(int i = 0; i < g_repeat; i++) { + CALL_SUBTEST( vectorProd(VectorXf(5)) ); + CALL_SUBTEST( vectorProd(VectorXd(10)) ); + CALL_SUBTEST( vectorProd(VectorXf(33)) ); + } +} |