1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2008 Benoit Jacob <jacob.benoit.1 (at) gmail.com> 5 // Copyright (C) 2010 Hauke Heibel <hauke.heibel (at) gmail.com> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #include "main.h" 12 13 #include <Eigen/StdList> 14 #include <Eigen/Geometry> 15 16 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Vector4f) 17 18 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Matrix2f) 19 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Matrix4f) 20 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Matrix4d) 21 22 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Affine3f) 23 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Affine3d) 24 25 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Quaternionf) 26 EIGEN_DEFINE_STL_LIST_SPECIALIZATION(Quaterniond) 27 28 template <class Container, class Position> 29 typename Container::iterator get(Container & c, Position position) 30 { 31 typename Container::iterator it = c.begin(); 32 std::advance(it, position); 33 return it; 34 } 35 36 template <class Container, class Position, class Value> 37 void set(Container & c, Position position, const Value & value) 38 { 39 typename Container::iterator it = c.begin(); 40 std::advance(it, position); 41 *it = value; 42 } 43 44 template<typename MatrixType> 45 void check_stdlist_matrix(const MatrixType& m) 46 { 47 typename MatrixType::Index rows = m.rows(); 48 typename MatrixType::Index cols = m.cols(); 49 MatrixType x = MatrixType::Random(rows,cols), y = MatrixType::Random(rows,cols); 50 std::list<MatrixType> v(10, MatrixType(rows,cols)), w(20, y); 51 typename std::list<MatrixType>::iterator itv = get(v, 5); 52 typename std::list<MatrixType>::iterator itw = get(w, 6); 53 *itv = x; 54 *itw = *itv; 55 VERIFY_IS_APPROX(*itw, *itv); 56 v = w; 57 itv = v.begin(); 58 itw = w.begin(); 59 for(int i = 0; i < 20; i++) 60 { 61 VERIFY_IS_APPROX(*itw, *itv); 62 ++itv; 63 ++itw; 64 } 65 66 v.resize(21); 67 set(v, 20, x); 68 VERIFY_IS_APPROX(*get(v, 20), x); 69 v.resize(22,y); 70 VERIFY_IS_APPROX(*get(v, 21), y); 71 v.push_back(x); 72 VERIFY_IS_APPROX(*get(v, 22), x); 73 74 // do a lot of push_back such that the list gets internally resized 75 // (with memory reallocation) 76 MatrixType* ref = &(*get(w, 0)); 77 for(int i=0; i<30 || ((ref==&(*get(w, 0))) && i<300); ++i) 78 v.push_back(*get(w, i%w.size())); 79 for(unsigned int i=23; i<v.size(); ++i) 80 { 81 VERIFY((*get(v, i))==(*get(w, (i-23)%w.size()))); 82 } 83 } 84 85 template<typename TransformType> 86 void check_stdlist_transform(const TransformType&) 87 { 88 typedef typename TransformType::MatrixType MatrixType; 89 TransformType x(MatrixType::Random()), y(MatrixType::Random()); 90 std::list<TransformType> v(10), w(20, y); 91 typename std::list<TransformType>::iterator itv = get(v, 5); 92 typename std::list<TransformType>::iterator itw = get(w, 6); 93 *itv = x; 94 *itw = *itv; 95 VERIFY_IS_APPROX(*itw, *itv); 96 v = w; 97 itv = v.begin(); 98 itw = w.begin(); 99 for(int i = 0; i < 20; i++) 100 { 101 VERIFY_IS_APPROX(*itw, *itv); 102 ++itv; 103 ++itw; 104 } 105 106 v.resize(21); 107 set(v, 20, x); 108 VERIFY_IS_APPROX(*get(v, 20), x); 109 v.resize(22,y); 110 VERIFY_IS_APPROX(*get(v, 21), y); 111 v.push_back(x); 112 VERIFY_IS_APPROX(*get(v, 22), x); 113 114 // do a lot of push_back such that the list gets internally resized 115 // (with memory reallocation) 116 TransformType* ref = &(*get(w, 0)); 117 for(int i=0; i<30 || ((ref==&(*get(w, 0))) && i<300); ++i) 118 v.push_back(*get(w, i%w.size())); 119 for(unsigned int i=23; i<v.size(); ++i) 120 { 121 VERIFY(get(v, i)->matrix()==get(w, (i-23)%w.size())->matrix()); 122 } 123 } 124 125 template<typename QuaternionType> 126 void check_stdlist_quaternion(const QuaternionType&) 127 { 128 typedef typename QuaternionType::Coefficients Coefficients; 129 QuaternionType x(Coefficients::Random()), y(Coefficients::Random()); 130 std::list<QuaternionType> v(10), w(20, y); 131 typename std::list<QuaternionType>::iterator itv = get(v, 5); 132 typename std::list<QuaternionType>::iterator itw = get(w, 6); 133 *itv = x; 134 *itw = *itv; 135 VERIFY_IS_APPROX(*itw, *itv); 136 v = w; 137 itv = v.begin(); 138 itw = w.begin(); 139 for(int i = 0; i < 20; i++) 140 { 141 VERIFY_IS_APPROX(*itw, *itv); 142 ++itv; 143 ++itw; 144 } 145 146 v.resize(21); 147 set(v, 20, x); 148 VERIFY_IS_APPROX(*get(v, 20), x); 149 v.resize(22,y); 150 VERIFY_IS_APPROX(*get(v, 21), y); 151 v.push_back(x); 152 VERIFY_IS_APPROX(*get(v, 22), x); 153 154 // do a lot of push_back such that the list gets internally resized 155 // (with memory reallocation) 156 QuaternionType* ref = &(*get(w, 0)); 157 for(int i=0; i<30 || ((ref==&(*get(w, 0))) && i<300); ++i) 158 v.push_back(*get(w, i%w.size())); 159 for(unsigned int i=23; i<v.size(); ++i) 160 { 161 VERIFY(get(v, i)->coeffs()==get(w, (i-23)%w.size())->coeffs()); 162 } 163 } 164 165 void test_stdlist_overload() 166 { 167 // some non vectorizable fixed sizes 168 CALL_SUBTEST_1(check_stdlist_matrix(Vector2f())); 169 CALL_SUBTEST_1(check_stdlist_matrix(Matrix3f())); 170 CALL_SUBTEST_2(check_stdlist_matrix(Matrix3d())); 171 172 // some vectorizable fixed sizes 173 CALL_SUBTEST_1(check_stdlist_matrix(Matrix2f())); 174 CALL_SUBTEST_1(check_stdlist_matrix(Vector4f())); 175 CALL_SUBTEST_1(check_stdlist_matrix(Matrix4f())); 176 CALL_SUBTEST_2(check_stdlist_matrix(Matrix4d())); 177 178 // some dynamic sizes 179 CALL_SUBTEST_3(check_stdlist_matrix(MatrixXd(1,1))); 180 CALL_SUBTEST_3(check_stdlist_matrix(VectorXd(20))); 181 CALL_SUBTEST_3(check_stdlist_matrix(RowVectorXf(20))); 182 CALL_SUBTEST_3(check_stdlist_matrix(MatrixXcf(10,10))); 183 184 // some Transform 185 CALL_SUBTEST_4(check_stdlist_transform(Affine2f())); // does not need the specialization (2+1)^2 = 9 186 CALL_SUBTEST_4(check_stdlist_transform(Affine3f())); 187 CALL_SUBTEST_4(check_stdlist_transform(Affine3d())); 188 189 // some Quaternion 190 CALL_SUBTEST_5(check_stdlist_quaternion(Quaternionf())); 191 CALL_SUBTEST_5(check_stdlist_quaternion(Quaterniond())); 192 } 193