diff options
author | Benoit Jacob <jacob.benoit.1@gmail.com> | 2011-03-06 20:59:25 -0500 |
---|---|---|
committer | Benoit Jacob <jacob.benoit.1@gmail.com> | 2011-03-06 20:59:25 -0500 |
commit | bfcad536e85cb189936eef3dd8e62597281ab32b (patch) | |
tree | a7e3d56cafbfeb652fda1f0d870967a097ac4398 /test/jacobisvd.cpp | |
parent | b464fc19bcaa50051fac345ea2f4ccdf52c3370e (diff) |
* bug #206: correctly forward computationOptions and work towards avoiding mallocs after preallocation, with unit test.
* added EIGEN_RUNTIME_NO_MALLOC and new set_is_malloc_allowed() function to implement that test
Diffstat (limited to 'test/jacobisvd.cpp')
-rw-r--r-- | test/jacobisvd.cpp | 46 |
1 files changed, 46 insertions, 0 deletions
diff --git a/test/jacobisvd.cpp b/test/jacobisvd.cpp index 72681bb73..907b290af 100644 --- a/test/jacobisvd.cpp +++ b/test/jacobisvd.cpp @@ -23,6 +23,9 @@ // License and a copy of the GNU General Public License along with // Eigen. If not, see <http://www.gnu.org/licenses/>. +// discard stack allocation as that too bypasses malloc +#define EIGEN_STACK_ALLOCATION_LIMIT 0 +#define EIGEN_RUNTIME_NO_MALLOC #include "main.h" #include <Eigen/SVD> @@ -241,6 +244,46 @@ void jacobisvd_inf_nan() svd.compute(m, ComputeFullU | ComputeFullV); } +void jacobisvd_preallocate() +{ + Vector3f v(3.f, 2.f, 1.f); + MatrixXf m = v.asDiagonal(); + + internal::set_is_malloc_allowed(false); + VERIFY_RAISES_ASSERT(VectorXf v(10);) + JacobiSVD<MatrixXf> svd; + internal::set_is_malloc_allowed(true); + svd.compute(m); + VERIFY_IS_APPROX(svd.singularValues(), v); + + JacobiSVD<MatrixXf> svd2(3,3); + internal::set_is_malloc_allowed(false); + svd2.compute(m); + internal::set_is_malloc_allowed(true); + VERIFY_IS_APPROX(svd2.singularValues(), v); + VERIFY_RAISES_ASSERT(svd2.matrixU()); + VERIFY_RAISES_ASSERT(svd2.matrixV()); + svd2.compute(m, ComputeFullU | ComputeFullV); + VERIFY_IS_APPROX(svd2.matrixU(), Matrix3f::Identity()); + VERIFY_IS_APPROX(svd2.matrixV(), Matrix3f::Identity()); + internal::set_is_malloc_allowed(false); + svd2.compute(m); + internal::set_is_malloc_allowed(true); + + JacobiSVD<MatrixXf> svd3(3,3,ComputeFullU|ComputeFullV); + internal::set_is_malloc_allowed(false); + svd2.compute(m); + internal::set_is_malloc_allowed(true); + VERIFY_IS_APPROX(svd2.singularValues(), v); + VERIFY_IS_APPROX(svd2.matrixU(), Matrix3f::Identity()); + VERIFY_IS_APPROX(svd2.matrixV(), Matrix3f::Identity()); + internal::set_is_malloc_allowed(false); + svd2.compute(m, ComputeFullU|ComputeFullV); + internal::set_is_malloc_allowed(true); + + +} + void test_jacobisvd() { CALL_SUBTEST_3(( jacobisvd_verify_assert(Matrix3f()) )); @@ -290,4 +333,7 @@ void test_jacobisvd() // Test problem size constructors CALL_SUBTEST_7( JacobiSVD<MatrixXf>(10,10) ); + + // Check that preallocation avoids subsequent mallocs + CALL_SUBTEST_9( jacobisvd_preallocate() ); } |