From 7f5256f628fbc860a4f85fae999f22e8cb7d2837 Mon Sep 17 00:00:00 2001 From: Hauke Heibel Date: Thu, 3 Sep 2009 17:27:51 +0200 Subject: [PATCH] Added conservativeResize + unit test. --- Eigen/src/Core/Matrix.h | 61 ++++++++++++++ Eigen/src/Core/util/StaticAssert.h | 6 ++ test/CMakeLists.txt | 1 + test/conservative_resize.cpp | 129 +++++++++++++++++++++++++++++ 4 files changed, 197 insertions(+) create mode 100644 test/conservative_resize.cpp diff --git a/Eigen/src/Core/Matrix.h b/Eigen/src/Core/Matrix.h index 2fc38c812..70862ad11 100644 --- a/Eigen/src/Core/Matrix.h +++ b/Eigen/src/Core/Matrix.h @@ -321,6 +321,67 @@ class Matrix else resize(other.rows(), other.cols()); } + /** Resizes \c *this to a \a rows x \a cols matrix while leaving old values of *this untouched. + * + * This method is intended for dynamic-size matrices, although it is legal to call it on any + * matrix as long as fixed dimensions are left unchanged. If you only want to change the number + * of rows and/or of columns, you can use conservativeResize(NoChange_t, int), + * conservativeResize(int, NoChange_t). + * + * The top-left part of the resized matrix will be the same as the overlapping top-left corner + * of *this. In case values need to be appended to the matrix they will be uninitialized per + * default and set to zero when init_with_zero is set to true. + */ + inline void conservativeResize(int rows, int cols, bool init_with_zero = false) + { + // Note: Here is space for improvement. Basically, for conservativeResize(int,int), + // neither RowsAtCompileTime or ColsAtCompileTime must be Dynamic. If only one of the + // dimensions is dynamic, one could use either conservativeResize(int rows, NoChange_t) or + // conservativeResize(NoChange_t, int cols). For these methods new static asserts like + // EIGEN_STATIC_ASSERT_DYNAMIC_ROWS and EIGEN_STATIC_ASSERT_DYNAMIC_COLS would be good. + EIGEN_STATIC_ASSERT_DYNAMIC_SIZE(Matrix) + PlainMatrixType tmp = init_with_zero ? PlainMatrixType::Zero(rows, cols) : PlainMatrixType(rows,cols); + const int common_rows = std::min(rows, this->rows()); + const int common_cols = std::min(cols, this->cols()); + tmp.block(0,0,common_rows,common_cols) = this->block(0,0,common_rows,common_cols); + this->derived().swap(tmp); + } + + EIGEN_STRONG_INLINE void conservativeResize(int rows, NoChange_t, bool init_with_zero = false) + { + // Note: see the comment in conservativeResize(int,int,bool) + conservativeResize(rows, cols(), init_with_zero); + } + + EIGEN_STRONG_INLINE void conservativeResize(NoChange_t, int cols, bool init_with_zero = false) + { + // Note: see the comment in conservativeResize(int,int,bool) + conservativeResize(rows(), cols, init_with_zero); + } + + /** Resizes \c *this to a vector of length \a size while retaining old values of *this. + * + * \only_for_vectors. This method does not work for + * partially dynamic matrices when the static dimension is anything other + * than 1. For example it will not work with Matrix. + * + * When values are appended, they will be uninitialized per default and set + * to zero when init_with_zero is set to true. + */ + inline void conservativeResize(int size, bool init_with_zero = false) + { + EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix) + EIGEN_STATIC_ASSERT_DYNAMIC_SIZE(Matrix) + + if (RowsAtCompileTime == 1 || ColsAtCompileTime == 1) + { + PlainMatrixType tmp = init_with_zero ? PlainMatrixType::Zero(size) : PlainMatrixType(size); + const int common_size = std::min(this->size(),size); + tmp.segment(0,common_size) = this->segment(0,common_size); + this->derived().swap(tmp); + } + } + /** Copies the value of the expression \a other into \c *this with automatic resizing. * * *this might be resized to match the dimensions of \a other. If *this was a null matrix (not already initialized), diff --git a/Eigen/src/Core/util/StaticAssert.h b/Eigen/src/Core/util/StaticAssert.h index e1ccd58e7..73d302fda 100644 --- a/Eigen/src/Core/util/StaticAssert.h +++ b/Eigen/src/Core/util/StaticAssert.h @@ -62,6 +62,7 @@ THIS_METHOD_IS_ONLY_FOR_MATRICES_OF_A_SPECIFIC_SIZE, YOU_MADE_A_PROGRAMMING_MISTAKE, YOU_CALLED_A_FIXED_SIZE_METHOD_ON_A_DYNAMIC_SIZE_MATRIX_OR_VECTOR, + YOU_CALLED_A_DYNAMIC_SIZE_METHOD_ON_A_FIXED_SIZE_MATRIX_OR_VECTOR, UNALIGNED_LOAD_AND_STORE_OPERATIONS_UNIMPLEMENTED_ON_ALTIVEC, NUMERIC_TYPE_MUST_BE_FLOATING_POINT, NUMERIC_TYPE_MUST_BE_REAL, @@ -114,6 +115,11 @@ EIGEN_STATIC_ASSERT(TYPE::SizeAtCompileTime!=Eigen::Dynamic, \ YOU_CALLED_A_FIXED_SIZE_METHOD_ON_A_DYNAMIC_SIZE_MATRIX_OR_VECTOR) +// static assertion failing if the type \a TYPE is not dynamic-size +#define EIGEN_STATIC_ASSERT_DYNAMIC_SIZE(TYPE) \ + EIGEN_STATIC_ASSERT(TYPE::SizeAtCompileTime==Eigen::Dynamic, \ + YOU_CALLED_A_DYNAMIC_SIZE_METHOD_ON_A_FIXED_SIZE_MATRIX_OR_VECTOR) + // static assertion failing if the type \a TYPE is not a vector type of the given size #define EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(TYPE, SIZE) \ EIGEN_STATIC_ASSERT(TYPE::IsVectorAtCompileTime && TYPE::SizeAtCompileTime==SIZE, \ diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 6bbf82144..03fbd48fc 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -149,6 +149,7 @@ ei_add_test(sparse_solvers " " "${SPARSE_LIBS}") ei_add_test(umeyama) ei_add_test(householder) ei_add_test(swap) +ei_add_test(conservative_resize) ei_add_property(EIGEN_TESTING_SUMMARY "CXX: ${CMAKE_CXX_COMPILER}\n") if(CMAKE_COMPILER_IS_GNUCXX) diff --git a/test/conservative_resize.cpp b/test/conservative_resize.cpp new file mode 100644 index 000000000..f0d025283 --- /dev/null +++ b/test/conservative_resize.cpp @@ -0,0 +1,129 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2009 Hauke Heibel +// +// 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 or1 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 . + +#include "main.h" + +#include +#include + +using namespace Eigen; + +template +void run_matrix_tests() +{ + typedef Matrix MatrixType; + + MatrixType m, n; + + // boundary cases ... + m = n = MatrixType::Random(50,50); + m.conservativeResize(1,50); + VERIFY_IS_APPROX(m, n.block(0,0,1,50)); + + m = n = MatrixType::Random(50,50); + m.conservativeResize(50,1); + VERIFY_IS_APPROX(m, n.block(0,0,50,1)); + + m = n = MatrixType::Random(50,50); + m.conservativeResize(50,50); + VERIFY_IS_APPROX(m, n.block(0,0,50,50)); + + // random shrinking ... + for (int i=0; i<25; ++i) + { + const int rows = ei_random(1,50); + const int cols = ei_random(1,50); + m = n = MatrixType::Random(50,50); + m.conservativeResize(rows,cols); + VERIFY_IS_APPROX(m, n.block(0,0,rows,cols)); + } + + // random growing with zeroing ... + for (int i=0; i<25; ++i) + { + const int rows = ei_random(50,75); + const int cols = ei_random(50,75); + m = n = MatrixType::Random(50,50); + m.conservativeResize(rows,cols,true); + VERIFY_IS_APPROX(m.block(0,0,n.rows(),n.cols()), n); + VERIFY( rows<=50 || m.block(50,0,rows-50,cols).sum() == Scalar(0) ); + VERIFY( cols<=50 || m.block(0,50,rows,cols-50).sum() == Scalar(0) ); + } +} + +template +void run_vector_tests() +{ + typedef Matrix MatrixType; + + MatrixType m, n; + + // boundary cases ... + m = n = MatrixType::Random(50); + m.conservativeResize(1); + VERIFY_IS_APPROX(m, n.segment(0,1)); + + m = n = MatrixType::Random(50); + m.conservativeResize(50); + VERIFY_IS_APPROX(m, n.segment(0,50)); + + // random shrinking ... + for (int i=0; i<50; ++i) + { + const int size = ei_random(1,50); + m = n = MatrixType::Random(50); + m.conservativeResize(size); + VERIFY_IS_APPROX(m, n.segment(0,size)); + } + + // random growing with zeroing ... + for (int i=0; i<50; ++i) + { + const int size = ei_random(50,100); + m = n = MatrixType::Random(50); + m.conservativeResize(size,true); + VERIFY_IS_APPROX(m.segment(0,50), n); + VERIFY( size<=50 || m.segment(50,size-50).sum() == Scalar(0) ); + } +} + +void test_conservative_resize() +{ + run_matrix_tests(); + run_matrix_tests(); + run_matrix_tests(); + run_matrix_tests(); + run_matrix_tests(); + run_matrix_tests(); + run_matrix_tests, Eigen::RowMajor>(); + run_matrix_tests, Eigen::ColMajor>(); + run_matrix_tests, Eigen::RowMajor>(); + run_matrix_tests, Eigen::ColMajor>(); + + run_vector_tests(); + run_vector_tests(); + run_vector_tests(); + run_vector_tests >(); + run_vector_tests >(); +}