1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 20013 Gael Guennebaud <gael.guennebaud@inria.fr>
5//
6// This Source Code Form is subject to the terms of the Mozilla
7// Public License v. 2.0. If a copy of the MPL was not distributed
8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10// This unit test cannot be easily written to work with EIGEN_DEFAULT_TO_ROW_MAJOR
11#ifdef EIGEN_DEFAULT_TO_ROW_MAJOR
12#undef EIGEN_DEFAULT_TO_ROW_MAJOR
13#endif
14
15static int nb_temporaries;
16
17inline void on_temporary_creation(int) {
18  // here's a great place to set a breakpoint when debugging failures in this test!
19  nb_temporaries++;
20}
21
22
23#define EIGEN_DENSE_STORAGE_CTOR_PLUGIN { on_temporary_creation(size); }
24
25#include "main.h"
26
27#define VERIFY_EVALUATION_COUNT(XPR,N) {\
28    nb_temporaries = 0; \
29    XPR; \
30    if(nb_temporaries!=N) std::cerr << "nb_temporaries == " << nb_temporaries << "\n"; \
31    VERIFY( (#XPR) && nb_temporaries==N ); \
32  }
33
34
35// test Ref.h
36
37template<typename MatrixType> void ref_matrix(const MatrixType& m)
38{
39  typedef typename MatrixType::Index Index;
40  typedef typename MatrixType::Scalar Scalar;
41  typedef typename MatrixType::RealScalar RealScalar;
42  typedef Matrix<Scalar,Dynamic,Dynamic,MatrixType::Options> DynMatrixType;
43  typedef Matrix<RealScalar,Dynamic,Dynamic,MatrixType::Options> RealDynMatrixType;
44
45  typedef Ref<MatrixType> RefMat;
46  typedef Ref<DynMatrixType> RefDynMat;
47  typedef Ref<const DynMatrixType> ConstRefDynMat;
48  typedef Ref<RealDynMatrixType , 0, Stride<Dynamic,Dynamic> > RefRealMatWithStride;
49
50  Index rows = m.rows(), cols = m.cols();
51
52  MatrixType  m1 = MatrixType::Random(rows, cols),
53              m2 = m1;
54
55  Index i = internal::random<Index>(0,rows-1);
56  Index j = internal::random<Index>(0,cols-1);
57  Index brows = internal::random<Index>(1,rows-i);
58  Index bcols = internal::random<Index>(1,cols-j);
59
60  RefMat rm0 = m1;
61  VERIFY_IS_EQUAL(rm0, m1);
62  RefDynMat rm1 = m1;
63  VERIFY_IS_EQUAL(rm1, m1);
64  RefDynMat rm2 = m1.block(i,j,brows,bcols);
65  VERIFY_IS_EQUAL(rm2, m1.block(i,j,brows,bcols));
66  rm2.setOnes();
67  m2.block(i,j,brows,bcols).setOnes();
68  VERIFY_IS_EQUAL(m1, m2);
69
70  m2.block(i,j,brows,bcols).setRandom();
71  rm2 = m2.block(i,j,brows,bcols);
72  VERIFY_IS_EQUAL(m1, m2);
73
74
75  ConstRefDynMat rm3 = m1.block(i,j,brows,bcols);
76  m1.block(i,j,brows,bcols) *= 2;
77  m2.block(i,j,brows,bcols) *= 2;
78  VERIFY_IS_EQUAL(rm3, m2.block(i,j,brows,bcols));
79  RefRealMatWithStride rm4 = m1.real();
80  VERIFY_IS_EQUAL(rm4, m2.real());
81  rm4.array() += 1;
82  m2.real().array() += 1;
83  VERIFY_IS_EQUAL(m1, m2);
84}
85
86template<typename VectorType> void ref_vector(const VectorType& m)
87{
88  typedef typename VectorType::Index Index;
89  typedef typename VectorType::Scalar Scalar;
90  typedef typename VectorType::RealScalar RealScalar;
91  typedef Matrix<Scalar,Dynamic,1,VectorType::Options> DynMatrixType;
92  typedef Matrix<Scalar,Dynamic,Dynamic,ColMajor> MatrixType;
93  typedef Matrix<RealScalar,Dynamic,1,VectorType::Options> RealDynMatrixType;
94
95  typedef Ref<VectorType> RefMat;
96  typedef Ref<DynMatrixType> RefDynMat;
97  typedef Ref<const DynMatrixType> ConstRefDynMat;
98  typedef Ref<RealDynMatrixType , 0, InnerStride<> > RefRealMatWithStride;
99  typedef Ref<DynMatrixType , 0, InnerStride<> > RefMatWithStride;
100
101  Index size = m.size();
102
103  VectorType  v1 = VectorType::Random(size),
104              v2 = v1;
105  MatrixType mat1 = MatrixType::Random(size,size),
106             mat2 = mat1,
107             mat3 = MatrixType::Random(size,size);
108
109  Index i = internal::random<Index>(0,size-1);
110  Index bsize = internal::random<Index>(1,size-i);
111
112  RefMat rm0 = v1;
113  VERIFY_IS_EQUAL(rm0, v1);
114  RefDynMat rv1 = v1;
115  VERIFY_IS_EQUAL(rv1, v1);
116  RefDynMat rv2 = v1.segment(i,bsize);
117  VERIFY_IS_EQUAL(rv2, v1.segment(i,bsize));
118  rv2.setOnes();
119  v2.segment(i,bsize).setOnes();
120  VERIFY_IS_EQUAL(v1, v2);
121
122  v2.segment(i,bsize).setRandom();
123  rv2 = v2.segment(i,bsize);
124  VERIFY_IS_EQUAL(v1, v2);
125
126  ConstRefDynMat rm3 = v1.segment(i,bsize);
127  v1.segment(i,bsize) *= 2;
128  v2.segment(i,bsize) *= 2;
129  VERIFY_IS_EQUAL(rm3, v2.segment(i,bsize));
130
131  RefRealMatWithStride rm4 = v1.real();
132  VERIFY_IS_EQUAL(rm4, v2.real());
133  rm4.array() += 1;
134  v2.real().array() += 1;
135  VERIFY_IS_EQUAL(v1, v2);
136
137  RefMatWithStride rm5 = mat1.row(i).transpose();
138  VERIFY_IS_EQUAL(rm5, mat1.row(i).transpose());
139  rm5.array() += 1;
140  mat2.row(i).array() += 1;
141  VERIFY_IS_EQUAL(mat1, mat2);
142  rm5.noalias() = rm4.transpose() * mat3;
143  mat2.row(i) = v2.real().transpose() * mat3;
144  VERIFY_IS_APPROX(mat1, mat2);
145}
146
147template<typename PlainObjectType> void check_const_correctness(const PlainObjectType&)
148{
149  // verify that ref-to-const don't have LvalueBit
150  typedef typename internal::add_const<PlainObjectType>::type ConstPlainObjectType;
151  VERIFY( !(internal::traits<Ref<ConstPlainObjectType> >::Flags & LvalueBit) );
152  VERIFY( !(internal::traits<Ref<ConstPlainObjectType, Aligned> >::Flags & LvalueBit) );
153  VERIFY( !(Ref<ConstPlainObjectType>::Flags & LvalueBit) );
154  VERIFY( !(Ref<ConstPlainObjectType, Aligned>::Flags & LvalueBit) );
155}
156
157template<typename B>
158EIGEN_DONT_INLINE void call_ref_1(Ref<VectorXf> a, const B &b) { VERIFY_IS_EQUAL(a,b); }
159template<typename B>
160EIGEN_DONT_INLINE void call_ref_2(const Ref<const VectorXf>& a, const B &b) { VERIFY_IS_EQUAL(a,b); }
161template<typename B>
162EIGEN_DONT_INLINE void call_ref_3(Ref<VectorXf,0,InnerStride<> > a, const B &b) { VERIFY_IS_EQUAL(a,b); }
163template<typename B>
164EIGEN_DONT_INLINE void call_ref_4(const Ref<const VectorXf,0,InnerStride<> >& a, const B &b) { VERIFY_IS_EQUAL(a,b); }
165template<typename B>
166EIGEN_DONT_INLINE void call_ref_5(Ref<MatrixXf,0,OuterStride<> > a, const B &b) { VERIFY_IS_EQUAL(a,b); }
167template<typename B>
168EIGEN_DONT_INLINE void call_ref_6(const Ref<const MatrixXf,0,OuterStride<> >& a, const B &b) { VERIFY_IS_EQUAL(a,b); }
169template<typename B>
170EIGEN_DONT_INLINE void call_ref_7(Ref<Matrix<float,Dynamic,3> > a, const B &b) { VERIFY_IS_EQUAL(a,b); }
171
172void call_ref()
173{
174  VectorXcf ca  = VectorXcf::Random(10);
175  VectorXf a    = VectorXf::Random(10);
176  RowVectorXf b = RowVectorXf::Random(10);
177  MatrixXf A    = MatrixXf::Random(10,10);
178  RowVector3f c = RowVector3f::Random();
179  const VectorXf& ac(a);
180  VectorBlock<VectorXf> ab(a,0,3);
181  const VectorBlock<VectorXf> abc(a,0,3);
182
183
184  VERIFY_EVALUATION_COUNT( call_ref_1(a,a), 0);
185  VERIFY_EVALUATION_COUNT( call_ref_1(b,b.transpose()), 0);
186//   call_ref_1(ac,a<c);           // does not compile because ac is const
187  VERIFY_EVALUATION_COUNT( call_ref_1(ab,ab), 0);
188  VERIFY_EVALUATION_COUNT( call_ref_1(a.head(4),a.head(4)), 0);
189  VERIFY_EVALUATION_COUNT( call_ref_1(abc,abc), 0);
190  VERIFY_EVALUATION_COUNT( call_ref_1(A.col(3),A.col(3)), 0);
191//   call_ref_1(A.row(3),A.row(3));    // does not compile because innerstride!=1
192  VERIFY_EVALUATION_COUNT( call_ref_3(A.row(3),A.row(3).transpose()), 0);
193  VERIFY_EVALUATION_COUNT( call_ref_4(A.row(3),A.row(3).transpose()), 0);
194//   call_ref_1(a+a, a+a);          // does not compile for obvious reason
195
196  MatrixXf tmp = A*A.col(1);
197  VERIFY_EVALUATION_COUNT( call_ref_2(A*A.col(1), tmp), 1);     // evaluated into a temp
198  VERIFY_EVALUATION_COUNT( call_ref_2(ac.head(5),ac.head(5)), 0);
199  VERIFY_EVALUATION_COUNT( call_ref_2(ac,ac), 0);
200  VERIFY_EVALUATION_COUNT( call_ref_2(a,a), 0);
201  VERIFY_EVALUATION_COUNT( call_ref_2(ab,ab), 0);
202  VERIFY_EVALUATION_COUNT( call_ref_2(a.head(4),a.head(4)), 0);
203  tmp = a+a;
204  VERIFY_EVALUATION_COUNT( call_ref_2(a+a,tmp), 1);            // evaluated into a temp
205  VERIFY_EVALUATION_COUNT( call_ref_2(ca.imag(),ca.imag()), 1);      // evaluated into a temp
206
207  VERIFY_EVALUATION_COUNT( call_ref_4(ac.head(5),ac.head(5)), 0);
208  tmp = a+a;
209  VERIFY_EVALUATION_COUNT( call_ref_4(a+a,tmp), 1);           // evaluated into a temp
210  VERIFY_EVALUATION_COUNT( call_ref_4(ca.imag(),ca.imag()), 0);
211
212  VERIFY_EVALUATION_COUNT( call_ref_5(a,a), 0);
213  VERIFY_EVALUATION_COUNT( call_ref_5(a.head(3),a.head(3)), 0);
214  VERIFY_EVALUATION_COUNT( call_ref_5(A,A), 0);
215//   call_ref_5(A.transpose(),A.transpose());   // does not compile because storage order does not match
216  VERIFY_EVALUATION_COUNT( call_ref_5(A.block(1,1,2,2),A.block(1,1,2,2)), 0);
217  VERIFY_EVALUATION_COUNT( call_ref_5(b,b), 0);             // storage order do not match, but this is a degenerate case that should work
218  VERIFY_EVALUATION_COUNT( call_ref_5(a.row(3),a.row(3)), 0);
219
220  VERIFY_EVALUATION_COUNT( call_ref_6(a,a), 0);
221  VERIFY_EVALUATION_COUNT( call_ref_6(a.head(3),a.head(3)), 0);
222  VERIFY_EVALUATION_COUNT( call_ref_6(A.row(3),A.row(3)), 1);           // evaluated into a temp thouth it could be avoided by viewing it as a 1xn matrix
223  tmp = A+A;
224  VERIFY_EVALUATION_COUNT( call_ref_6(A+A,tmp), 1);                // evaluated into a temp
225  VERIFY_EVALUATION_COUNT( call_ref_6(A,A), 0);
226  VERIFY_EVALUATION_COUNT( call_ref_6(A.transpose(),A.transpose()), 1);      // evaluated into a temp because the storage orders do not match
227  VERIFY_EVALUATION_COUNT( call_ref_6(A.block(1,1,2,2),A.block(1,1,2,2)), 0);
228
229  VERIFY_EVALUATION_COUNT( call_ref_7(c,c), 0);
230}
231
232typedef Matrix<double,Dynamic,Dynamic,RowMajor> RowMatrixXd;
233int test_ref_overload_fun1(Ref<MatrixXd> )       { return 1; }
234int test_ref_overload_fun1(Ref<RowMatrixXd> )    { return 2; }
235int test_ref_overload_fun1(Ref<MatrixXf> )       { return 3; }
236
237int test_ref_overload_fun2(Ref<const MatrixXd> ) { return 4; }
238int test_ref_overload_fun2(Ref<const MatrixXf> ) { return 5; }
239
240// See also bug 969
241void test_ref_overloads()
242{
243  MatrixXd Ad, Bd;
244  RowMatrixXd rAd, rBd;
245  VERIFY( test_ref_overload_fun1(Ad)==1 );
246  VERIFY( test_ref_overload_fun1(rAd)==2 );
247
248  MatrixXf Af, Bf;
249  VERIFY( test_ref_overload_fun2(Ad)==4 );
250  VERIFY( test_ref_overload_fun2(Ad+Bd)==4 );
251  VERIFY( test_ref_overload_fun2(Af+Bf)==5 );
252}
253
254void test_ref()
255{
256  for(int i = 0; i < g_repeat; i++) {
257    CALL_SUBTEST_1( ref_vector(Matrix<float, 1, 1>()) );
258    CALL_SUBTEST_1( check_const_correctness(Matrix<float, 1, 1>()) );
259    CALL_SUBTEST_2( ref_vector(Vector4d()) );
260    CALL_SUBTEST_2( check_const_correctness(Matrix4d()) );
261    CALL_SUBTEST_3( ref_vector(Vector4cf()) );
262    CALL_SUBTEST_4( ref_vector(VectorXcf(8)) );
263    CALL_SUBTEST_5( ref_vector(VectorXi(12)) );
264    CALL_SUBTEST_5( check_const_correctness(VectorXi(12)) );
265
266    CALL_SUBTEST_1( ref_matrix(Matrix<float, 1, 1>()) );
267    CALL_SUBTEST_2( ref_matrix(Matrix4d()) );
268    CALL_SUBTEST_1( ref_matrix(Matrix<float,3,5>()) );
269    CALL_SUBTEST_4( ref_matrix(MatrixXcf(internal::random<int>(1,10),internal::random<int>(1,10))) );
270    CALL_SUBTEST_4( ref_matrix(Matrix<std::complex<double>,10,15>()) );
271    CALL_SUBTEST_5( ref_matrix(MatrixXi(internal::random<int>(1,10),internal::random<int>(1,10))) );
272    CALL_SUBTEST_6( call_ref() );
273  }
274
275  CALL_SUBTEST_7( test_ref_overloads() );
276}
277