1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
5 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@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 #define EIGEN2_SUPPORT
12 #define EIGEN_NO_EIGEN2_DEPRECATED_WARNING
13
14 #define EIGEN_NO_STATIC_ASSERT
15 #include "main.h"
16 #include <functional>
17
18 #ifdef min
19 #undef min
20 #endif
21
22 #ifdef max
23 #undef max
24 #endif
25
26 using namespace std;
27
28 template<typename Scalar> struct AddIfNull {
operator ()AddIfNull29 const Scalar operator() (const Scalar a, const Scalar b) const {return a<=1e-3 ? b : a;}
30 enum { Cost = NumTraits<Scalar>::AddCost };
31 };
32
33 template<typename MatrixType>
34 typename Eigen::internal::enable_if<!NumTraits<typename MatrixType::Scalar>::IsInteger,typename MatrixType::Scalar>::type
cwiseops_real_only(MatrixType & m1,MatrixType & m2,MatrixType & m3,MatrixType & mones)35 cwiseops_real_only(MatrixType& m1, MatrixType& m2, MatrixType& m3, MatrixType& mones)
36 {
37 typedef typename MatrixType::Scalar Scalar;
38 typedef typename NumTraits<Scalar>::Real RealScalar;
39
40 VERIFY_IS_APPROX(m1.cwise() / m2, m1.cwise() * (m2.cwise().inverse()));
41 m3 = m1.cwise().abs().cwise().sqrt();
42 VERIFY_IS_APPROX(m3.cwise().square(), m1.cwise().abs());
43 VERIFY_IS_APPROX(m1.cwise().square().cwise().sqrt(), m1.cwise().abs());
44 VERIFY_IS_APPROX(m1.cwise().abs().cwise().log().cwise().exp() , m1.cwise().abs());
45
46 VERIFY_IS_APPROX(m1.cwise().pow(2), m1.cwise().square());
47 m3 = (m1.cwise().abs().cwise()<=RealScalar(0.01)).select(mones,m1);
48 VERIFY_IS_APPROX(m3.cwise().pow(-1), m3.cwise().inverse());
49 m3 = m1.cwise().abs();
50 VERIFY_IS_APPROX(m3.cwise().pow(RealScalar(0.5)), m3.cwise().sqrt());
51
52 // VERIFY_IS_APPROX(m1.cwise().tan(), m1.cwise().sin().cwise() / m1.cwise().cos());
53 VERIFY_IS_APPROX(mones, m1.cwise().sin().cwise().square() + m1.cwise().cos().cwise().square());
54 m3 = m1;
55 m3.cwise() /= m2;
56 VERIFY_IS_APPROX(m3, m1.cwise() / m2);
57
58 return Scalar(0);
59 }
60
61 template<typename MatrixType>
62 typename Eigen::internal::enable_if<NumTraits<typename MatrixType::Scalar>::IsInteger,typename MatrixType::Scalar>::type
cwiseops_real_only(MatrixType &,MatrixType &,MatrixType &,MatrixType &)63 cwiseops_real_only(MatrixType& , MatrixType& , MatrixType& , MatrixType& )
64 {
65 return 0;
66 }
67
cwiseops(const MatrixType & m)68 template<typename MatrixType> void cwiseops(const MatrixType& m)
69 {
70 typedef typename MatrixType::Index Index;
71 typedef typename MatrixType::Scalar Scalar;
72 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType;
73
74 Index rows = m.rows();
75 Index cols = m.cols();
76
77 MatrixType m1 = MatrixType::Random(rows, cols),
78 m1bis = m1,
79 m2 = MatrixType::Random(rows, cols),
80 m3(rows, cols),
81 m4(rows, cols),
82 mzero = MatrixType::Zero(rows, cols),
83 mones = MatrixType::Ones(rows, cols),
84 identity = Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime>
85 ::Identity(rows, rows);
86 VectorType vzero = VectorType::Zero(rows),
87 vones = VectorType::Ones(rows),
88 v3(rows);
89
90 Index r = internal::random<Index>(0, rows-1),
91 c = internal::random<Index>(0, cols-1);
92
93 Scalar s1 = internal::random<Scalar>();
94
95 // test Zero, Ones, Constant, and the set* variants
96 m3 = MatrixType::Constant(rows, cols, s1);
97 for (int j=0; j<cols; ++j)
98 for (int i=0; i<rows; ++i)
99 {
100 VERIFY_IS_APPROX(mzero(i,j), Scalar(0));
101 VERIFY_IS_APPROX(mones(i,j), Scalar(1));
102 VERIFY_IS_APPROX(m3(i,j), s1);
103 }
104 VERIFY(mzero.isZero());
105 VERIFY(mones.isOnes());
106 VERIFY(m3.isConstant(s1));
107 VERIFY(identity.isIdentity());
108 VERIFY_IS_APPROX(m4.setConstant(s1), m3);
109 VERIFY_IS_APPROX(m4.setConstant(rows,cols,s1), m3);
110 VERIFY_IS_APPROX(m4.setZero(), mzero);
111 VERIFY_IS_APPROX(m4.setZero(rows,cols), mzero);
112 VERIFY_IS_APPROX(m4.setOnes(), mones);
113 VERIFY_IS_APPROX(m4.setOnes(rows,cols), mones);
114 m4.fill(s1);
115 VERIFY_IS_APPROX(m4, m3);
116
117 VERIFY_IS_APPROX(v3.setConstant(rows, s1), VectorType::Constant(rows,s1));
118 VERIFY_IS_APPROX(v3.setZero(rows), vzero);
119 VERIFY_IS_APPROX(v3.setOnes(rows), vones);
120
121 m2 = m2.template binaryExpr<AddIfNull<Scalar> >(mones);
122
123 VERIFY_IS_APPROX(m1.cwise().pow(2), m1.cwise().abs2());
124 VERIFY_IS_APPROX(m1.cwise().pow(2), m1.cwise().square());
125 VERIFY_IS_APPROX(m1.cwise().pow(3), m1.cwise().cube());
126
127 VERIFY_IS_APPROX(m1 + mones, m1.cwise()+Scalar(1));
128 VERIFY_IS_APPROX(m1 - mones, m1.cwise()-Scalar(1));
129 m3 = m1; m3.cwise() += 1;
130 VERIFY_IS_APPROX(m1 + mones, m3);
131 m3 = m1; m3.cwise() -= 1;
132 VERIFY_IS_APPROX(m1 - mones, m3);
133
134 VERIFY_IS_APPROX(m2, m2.cwise() * mones);
135 VERIFY_IS_APPROX(m1.cwise() * m2, m2.cwise() * m1);
136 m3 = m1;
137 m3.cwise() *= m2;
138 VERIFY_IS_APPROX(m3, m1.cwise() * m2);
139
140 VERIFY_IS_APPROX(mones, m2.cwise()/m2);
141
142 // check min
143 VERIFY_IS_APPROX( m1.cwise().min(m2), m2.cwise().min(m1) );
144 VERIFY_IS_APPROX( m1.cwise().min(m1+mones), m1 );
145 VERIFY_IS_APPROX( m1.cwise().min(m1-mones), m1-mones );
146
147 // check max
148 VERIFY_IS_APPROX( m1.cwise().max(m2), m2.cwise().max(m1) );
149 VERIFY_IS_APPROX( m1.cwise().max(m1-mones), m1 );
150 VERIFY_IS_APPROX( m1.cwise().max(m1+mones), m1+mones );
151
152 VERIFY( (m1.cwise() == m1).all() );
153 VERIFY( (m1.cwise() != m2).any() );
154 VERIFY(!(m1.cwise() == (m1+mones)).any() );
155 if (rows*cols>1)
156 {
157 m3 = m1;
158 m3(r,c) += 1;
159 VERIFY( (m1.cwise() == m3).any() );
160 VERIFY( !(m1.cwise() == m3).all() );
161 }
162 VERIFY( (m1.cwise().min(m2).cwise() <= m2).all() );
163 VERIFY( (m1.cwise().max(m2).cwise() >= m2).all() );
164 VERIFY( (m1.cwise().min(m2).cwise() < (m1+mones)).all() );
165 VERIFY( (m1.cwise().max(m2).cwise() > (m1-mones)).all() );
166
167 VERIFY( (m1.cwise()<m1.unaryExpr(bind2nd(plus<Scalar>(), Scalar(1)))).all() );
168 VERIFY( !(m1.cwise()<m1bis.unaryExpr(bind2nd(minus<Scalar>(), Scalar(1)))).all() );
169 VERIFY( !(m1.cwise()>m1bis.unaryExpr(bind2nd(plus<Scalar>(), Scalar(1)))).any() );
170
171 cwiseops_real_only(m1, m2, m3, mones);
172 }
173
test_cwiseop()174 void test_cwiseop()
175 {
176 for(int i = 0; i < g_repeat ; i++) {
177 CALL_SUBTEST_1( cwiseops(Matrix<float, 1, 1>()) );
178 CALL_SUBTEST_2( cwiseops(Matrix4d()) );
179 CALL_SUBTEST_3( cwiseops(MatrixXf(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
180 CALL_SUBTEST_4( cwiseops(MatrixXf(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
181 CALL_SUBTEST_5( cwiseops(MatrixXi(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
182 CALL_SUBTEST_6( cwiseops(MatrixXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
183 }
184 }
185