• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2011 Benoit Jacob <jacob.benoit.1@gmail.com>
5 // Copyright (C) 2015 Gael Guennebaud <gael.guennebaud@inria.fr>
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 #define TEST_ENABLE_TEMPORARY_TRACKING
12 #define EIGEN_NO_STATIC_ASSERT
13 
14 #include "main.h"
15 
vectorwiseop_array(const ArrayType & m)16 template<typename ArrayType> void vectorwiseop_array(const ArrayType& m)
17 {
18   typedef typename ArrayType::Index Index;
19   typedef typename ArrayType::Scalar Scalar;
20   typedef Array<Scalar, ArrayType::RowsAtCompileTime, 1> ColVectorType;
21   typedef Array<Scalar, 1, ArrayType::ColsAtCompileTime> RowVectorType;
22 
23   Index rows = m.rows();
24   Index cols = m.cols();
25   Index r = internal::random<Index>(0, rows-1),
26         c = internal::random<Index>(0, cols-1);
27 
28   ArrayType m1 = ArrayType::Random(rows, cols),
29             m2(rows, cols),
30             m3(rows, cols);
31 
32   ColVectorType colvec = ColVectorType::Random(rows);
33   RowVectorType rowvec = RowVectorType::Random(cols);
34 
35   // test addition
36 
37   m2 = m1;
38   m2.colwise() += colvec;
39   VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
40   VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
41 
42   VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
43   VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
44 
45   m2 = m1;
46   m2.rowwise() += rowvec;
47   VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
48   VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
49 
50   VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
51   VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
52 
53   // test substraction
54 
55   m2 = m1;
56   m2.colwise() -= colvec;
57   VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
58   VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
59 
60   VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
61   VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
62 
63   m2 = m1;
64   m2.rowwise() -= rowvec;
65   VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
66   VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
67 
68   VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
69   VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
70 
71   // test multiplication
72 
73   m2 = m1;
74   m2.colwise() *= colvec;
75   VERIFY_IS_APPROX(m2, m1.colwise() * colvec);
76   VERIFY_IS_APPROX(m2.col(c), m1.col(c) * colvec);
77 
78   VERIFY_RAISES_ASSERT(m2.colwise() *= colvec.transpose());
79   VERIFY_RAISES_ASSERT(m1.colwise() * colvec.transpose());
80 
81   m2 = m1;
82   m2.rowwise() *= rowvec;
83   VERIFY_IS_APPROX(m2, m1.rowwise() * rowvec);
84   VERIFY_IS_APPROX(m2.row(r), m1.row(r) * rowvec);
85 
86   VERIFY_RAISES_ASSERT(m2.rowwise() *= rowvec.transpose());
87   VERIFY_RAISES_ASSERT(m1.rowwise() * rowvec.transpose());
88 
89   // test quotient
90 
91   m2 = m1;
92   m2.colwise() /= colvec;
93   VERIFY_IS_APPROX(m2, m1.colwise() / colvec);
94   VERIFY_IS_APPROX(m2.col(c), m1.col(c) / colvec);
95 
96   VERIFY_RAISES_ASSERT(m2.colwise() /= colvec.transpose());
97   VERIFY_RAISES_ASSERT(m1.colwise() / colvec.transpose());
98 
99   m2 = m1;
100   m2.rowwise() /= rowvec;
101   VERIFY_IS_APPROX(m2, m1.rowwise() / rowvec);
102   VERIFY_IS_APPROX(m2.row(r), m1.row(r) / rowvec);
103 
104   VERIFY_RAISES_ASSERT(m2.rowwise() /= rowvec.transpose());
105   VERIFY_RAISES_ASSERT(m1.rowwise() / rowvec.transpose());
106 
107   m2 = m1;
108   // yes, there might be an aliasing issue there but ".rowwise() /="
109   // is supposed to evaluate " m2.colwise().sum()" into a temporary to avoid
110   // evaluating the reduction multiple times
111   if(ArrayType::RowsAtCompileTime>2 || ArrayType::RowsAtCompileTime==Dynamic)
112   {
113     m2.rowwise() /= m2.colwise().sum();
114     VERIFY_IS_APPROX(m2, m1.rowwise() / m1.colwise().sum());
115   }
116 
117   // all/any
118   Array<bool,Dynamic,Dynamic> mb(rows,cols);
119   mb = (m1.real()<=0.7).colwise().all();
120   VERIFY( (mb.col(c) == (m1.real().col(c)<=0.7).all()).all() );
121   mb = (m1.real()<=0.7).rowwise().all();
122   VERIFY( (mb.row(r) == (m1.real().row(r)<=0.7).all()).all() );
123 
124   mb = (m1.real()>=0.7).colwise().any();
125   VERIFY( (mb.col(c) == (m1.real().col(c)>=0.7).any()).all() );
126   mb = (m1.real()>=0.7).rowwise().any();
127   VERIFY( (mb.row(r) == (m1.real().row(r)>=0.7).any()).all() );
128 }
129 
vectorwiseop_matrix(const MatrixType & m)130 template<typename MatrixType> void vectorwiseop_matrix(const MatrixType& m)
131 {
132   typedef typename MatrixType::Index Index;
133   typedef typename MatrixType::Scalar Scalar;
134   typedef typename NumTraits<Scalar>::Real RealScalar;
135   typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType;
136   typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType;
137   typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, 1> RealColVectorType;
138   typedef Matrix<RealScalar, 1, MatrixType::ColsAtCompileTime> RealRowVectorType;
139 
140   Index rows = m.rows();
141   Index cols = m.cols();
142   Index r = internal::random<Index>(0, rows-1),
143         c = internal::random<Index>(0, cols-1);
144 
145   MatrixType m1 = MatrixType::Random(rows, cols),
146             m2(rows, cols),
147             m3(rows, cols);
148 
149   ColVectorType colvec = ColVectorType::Random(rows);
150   RowVectorType rowvec = RowVectorType::Random(cols);
151   RealColVectorType rcres;
152   RealRowVectorType rrres;
153 
154   // test addition
155 
156   m2 = m1;
157   m2.colwise() += colvec;
158   VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
159   VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
160 
161   if(rows>1)
162   {
163     VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
164     VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
165   }
166 
167   m2 = m1;
168   m2.rowwise() += rowvec;
169   VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
170   VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
171 
172   if(cols>1)
173   {
174     VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
175     VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
176   }
177 
178   // test substraction
179 
180   m2 = m1;
181   m2.colwise() -= colvec;
182   VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
183   VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
184 
185   if(rows>1)
186   {
187     VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
188     VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
189   }
190 
191   m2 = m1;
192   m2.rowwise() -= rowvec;
193   VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
194   VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
195 
196   if(cols>1)
197   {
198     VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
199     VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
200   }
201 
202   // test norm
203   rrres = m1.colwise().norm();
204   VERIFY_IS_APPROX(rrres(c), m1.col(c).norm());
205   rcres = m1.rowwise().norm();
206   VERIFY_IS_APPROX(rcres(r), m1.row(r).norm());
207 
208   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().sum(), m1.colwise().template lpNorm<1>());
209   VERIFY_IS_APPROX(m1.cwiseAbs().rowwise().sum(), m1.rowwise().template lpNorm<1>());
210   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().maxCoeff(), m1.colwise().template lpNorm<Infinity>());
211   VERIFY_IS_APPROX(m1.cwiseAbs().rowwise().maxCoeff(), m1.rowwise().template lpNorm<Infinity>());
212 
213   // regression for bug 1158
214   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().sum().x(), m1.col(0).cwiseAbs().sum());
215 
216   // test normalized
217   m2 = m1.colwise().normalized();
218   VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
219   m2 = m1.rowwise().normalized();
220   VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
221 
222   // test normalize
223   m2 = m1;
224   m2.colwise().normalize();
225   VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
226   m2 = m1;
227   m2.rowwise().normalize();
228   VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
229 
230   // test with partial reduction of products
231   Matrix<Scalar,MatrixType::RowsAtCompileTime,MatrixType::RowsAtCompileTime> m1m1 = m1 * m1.transpose();
232   VERIFY_IS_APPROX( (m1 * m1.transpose()).colwise().sum(), m1m1.colwise().sum());
233   Matrix<Scalar,1,MatrixType::RowsAtCompileTime> tmp(rows);
234   VERIFY_EVALUATION_COUNT( tmp = (m1 * m1.transpose()).colwise().sum(), 1);
235 
236   m2 = m1.rowwise() - (m1.colwise().sum()/RealScalar(m1.rows())).eval();
237   m1 = m1.rowwise() - (m1.colwise().sum()/RealScalar(m1.rows()));
238   VERIFY_IS_APPROX( m1, m2 );
239   VERIFY_EVALUATION_COUNT( m2 = (m1.rowwise() - m1.colwise().sum()/RealScalar(m1.rows())), (MatrixType::RowsAtCompileTime!=1 ? 1 : 0) );
240 }
241 
test_vectorwiseop()242 void test_vectorwiseop()
243 {
244   CALL_SUBTEST_1( vectorwiseop_array(Array22cd()) );
245   CALL_SUBTEST_2( vectorwiseop_array(Array<double, 3, 2>()) );
246   CALL_SUBTEST_3( vectorwiseop_array(ArrayXXf(3, 4)) );
247   CALL_SUBTEST_4( vectorwiseop_matrix(Matrix4cf()) );
248   CALL_SUBTEST_5( vectorwiseop_matrix(Matrix<float,4,5>()) );
249   CALL_SUBTEST_6( vectorwiseop_matrix(MatrixXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
250   CALL_SUBTEST_7( vectorwiseop_matrix(VectorXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
251   CALL_SUBTEST_7( vectorwiseop_matrix(RowVectorXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
252 }
253