linearstructure.cpp
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
5 // Copyright (C) 2014 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 static bool g_called;
12 #define EIGEN_SCALAR_BINARY_OP_PLUGIN { g_called |= (!internal::is_same<LhsScalar,RhsScalar>::value); }
13 
14 #include "main.h"
15 
16 template<typename MatrixType> void linearStructure(const MatrixType& m)
17 {
18  using std::abs;
19  /* this test covers the following files:
20  CwiseUnaryOp.h, CwiseBinaryOp.h, SelfCwiseBinaryOp.h
21  */
22  typedef typename MatrixType::Scalar Scalar;
23  typedef typename MatrixType::RealScalar RealScalar;
24 
25  Index rows = m.rows();
26  Index cols = m.cols();
27 
28  // this test relies a lot on Random.h, and there's not much more that we can do
29  // to test it, hence I consider that we will have tested Random.h
30  MatrixType m1 = MatrixType::Random(rows, cols),
31  m2 = MatrixType::Random(rows, cols),
32  m3(rows, cols);
33 
34  Scalar s1 = internal::random<Scalar>();
35  while (abs(s1)<RealScalar(1e-3)) s1 = internal::random<Scalar>();
36 
37  Index r = internal::random<Index>(0, rows-1),
38  c = internal::random<Index>(0, cols-1);
39 
40  VERIFY_IS_APPROX(-(-m1), m1);
41  VERIFY_IS_APPROX(m1+m1, 2*m1);
42  VERIFY_IS_APPROX(m1+m2-m1, m2);
43  VERIFY_IS_APPROX(-m2+m1+m2, m1);
44  VERIFY_IS_APPROX(m1*s1, s1*m1);
45  VERIFY_IS_APPROX((m1+m2)*s1, s1*m1+s1*m2);
46  VERIFY_IS_APPROX((-m1+m2)*s1, -s1*m1+s1*m2);
47  m3 = m2; m3 += m1;
48  VERIFY_IS_APPROX(m3, m1+m2);
49  m3 = m2; m3 -= m1;
50  VERIFY_IS_APPROX(m3, m2-m1);
51  m3 = m2; m3 *= s1;
52  VERIFY_IS_APPROX(m3, s1*m2);
53  if(!NumTraits<Scalar>::IsInteger)
54  {
55  m3 = m2; m3 /= s1;
56  VERIFY_IS_APPROX(m3, m2/s1);
57  }
58 
59  // again, test operator() to check const-qualification
60  VERIFY_IS_APPROX((-m1)(r,c), -(m1(r,c)));
61  VERIFY_IS_APPROX((m1-m2)(r,c), (m1(r,c))-(m2(r,c)));
62  VERIFY_IS_APPROX((m1+m2)(r,c), (m1(r,c))+(m2(r,c)));
63  VERIFY_IS_APPROX((s1*m1)(r,c), s1*(m1(r,c)));
64  VERIFY_IS_APPROX((m1*s1)(r,c), (m1(r,c))*s1);
65  if(!NumTraits<Scalar>::IsInteger)
66  VERIFY_IS_APPROX((m1/s1)(r,c), (m1(r,c))/s1);
67 
68  // use .block to disable vectorization and compare to the vectorized version
69  VERIFY_IS_APPROX(m1+m1.block(0,0,rows,cols), m1+m1);
70  VERIFY_IS_APPROX(m1.cwiseProduct(m1.block(0,0,rows,cols)), m1.cwiseProduct(m1));
71  VERIFY_IS_APPROX(m1 - m1.block(0,0,rows,cols), m1 - m1);
72  VERIFY_IS_APPROX(m1.block(0,0,rows,cols) * s1, m1 * s1);
73 }
74 
75 // Make sure that complex * real and real * complex are properly optimized
76 template<typename MatrixType> void real_complex(DenseIndex rows = MatrixType::RowsAtCompileTime, DenseIndex cols = MatrixType::ColsAtCompileTime)
77 {
78  typedef typename MatrixType::Scalar Scalar;
79  typedef typename MatrixType::RealScalar RealScalar;
80 
81  RealScalar s = internal::random<RealScalar>();
82  MatrixType m1 = MatrixType::Random(rows, cols);
83 
84  g_called = false;
85  VERIFY_IS_APPROX(s*m1, Scalar(s)*m1);
86  VERIFY(g_called && "real * matrix<complex> not properly optimized");
87 
88  g_called = false;
89  VERIFY_IS_APPROX(m1*s, m1*Scalar(s));
90  VERIFY(g_called && "matrix<complex> * real not properly optimized");
91 
92  g_called = false;
93  VERIFY_IS_APPROX(m1/s, m1/Scalar(s));
94  VERIFY(g_called && "matrix<complex> / real not properly optimized");
95 
96  g_called = false;
97  VERIFY_IS_APPROX(s+m1.array(), Scalar(s)+m1.array());
98  VERIFY(g_called && "real + matrix<complex> not properly optimized");
99 
100  g_called = false;
101  VERIFY_IS_APPROX(m1.array()+s, m1.array()+Scalar(s));
102  VERIFY(g_called && "matrix<complex> + real not properly optimized");
103 
104  g_called = false;
105  VERIFY_IS_APPROX(s-m1.array(), Scalar(s)-m1.array());
106  VERIFY(g_called && "real - matrix<complex> not properly optimized");
107 
108  g_called = false;
109  VERIFY_IS_APPROX(m1.array()-s, m1.array()-Scalar(s));
110  VERIFY(g_called && "matrix<complex> - real not properly optimized");
111 }
112 
113 template<int>
114 void linearstructure_overflow()
115 {
116  // make sure that /=scalar and /scalar do not overflow
117  // rational: 1.0/4.94e-320 overflow, but m/4.94e-320 should not
118  Matrix4d m2, m3;
119  m3 = m2 = Matrix4d::Random()*1e-20;
120  m2 = m2 / 4.9e-320;
121  VERIFY_IS_APPROX(m2.cwiseQuotient(m2), Matrix4d::Ones());
122  m3 /= 4.9e-320;
123  VERIFY_IS_APPROX(m3.cwiseQuotient(m3), Matrix4d::Ones());
124 }
125 
126 EIGEN_DECLARE_TEST(linearstructure)
127 {
128  g_called = true;
129  VERIFY(g_called); // avoid `unneeded-internal-declaration` warning.
130  for(int i = 0; i < g_repeat; i++) {
131  CALL_SUBTEST_1( linearStructure(Matrix<float, 1, 1>()) );
132  CALL_SUBTEST_2( linearStructure(Matrix2f()) );
133  CALL_SUBTEST_3( linearStructure(Vector3d()) );
134  CALL_SUBTEST_4( linearStructure(Matrix4d()) );
135  CALL_SUBTEST_5( linearStructure(MatrixXcf(internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2), internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2))) );
136  CALL_SUBTEST_6( linearStructure(MatrixXf (internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
137  CALL_SUBTEST_7( linearStructure(MatrixXi (internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
138  CALL_SUBTEST_8( linearStructure(MatrixXcd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2), internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2))) );
139  CALL_SUBTEST_9( linearStructure(ArrayXXf (internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
140  CALL_SUBTEST_10( linearStructure(ArrayXXcf (internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
141 
142  CALL_SUBTEST_11( real_complex<Matrix4cd>() );
143  CALL_SUBTEST_11( real_complex<MatrixXcf>(10,10) );
144  CALL_SUBTEST_11( real_complex<ArrayXXcf>(10,10) );
145  }
146  CALL_SUBTEST_4( linearstructure_overflow<0>() );
147 }
m
Matrix3f m
Definition: AngleAxis_mimic_euler.cpp:1
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:46
CALL_SUBTEST_9
#define CALL_SUBTEST_9(FUNC)
Definition: split_test_helper.h:52
CALL_SUBTEST_6
#define CALL_SUBTEST_6(FUNC)
Definition: split_test_helper.h:34
CALL_SUBTEST_4
#define CALL_SUBTEST_4(FUNC)
Definition: split_test_helper.h:22
cols
int cols
Definition: Tutorial_commainit_02.cpp:1
CALL_SUBTEST_3
#define CALL_SUBTEST_3(FUNC)
Definition: split_test_helper.h:16
m2
MatrixType m2(n_dims)
CALL_SUBTEST_7
#define CALL_SUBTEST_7(FUNC)
Definition: split_test_helper.h:40
c
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
CALL_SUBTEST_11
#define CALL_SUBTEST_11(FUNC)
Definition: split_test_helper.h:64
MatrixType
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
Eigen::NumTraits
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
CALL_SUBTEST_10
#define CALL_SUBTEST_10(FUNC)
Definition: split_test_helper.h:58
VERIFY_IS_APPROX
#define VERIFY_IS_APPROX(a, b)
Definition: integer_types.cpp:15
CALL_SUBTEST_1
#define CALL_SUBTEST_1(FUNC)
Definition: split_test_helper.h:4
m1
Matrix3d m1
Definition: IOFormat.cpp:2
Eigen::g_repeat
static int g_repeat
Definition: main.h:169
Eigen::Index
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
CALL_SUBTEST_8
#define CALL_SUBTEST_8(FUNC)
Definition: split_test_helper.h:46
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
s
RealScalar s
Definition: level1_cplx_impl.h:126
RealScalar
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
g_called
static bool g_called
Definition: linearstructure.cpp:11
Eigen::DenseIndex
EIGEN_DEFAULT_DENSE_INDEX_TYPE DenseIndex
Definition: Meta.h:66
linearStructure
void linearStructure(const MatrixType &m)
Definition: linearstructure.cpp:16
CALL_SUBTEST_5
#define CALL_SUBTEST_5(FUNC)
Definition: split_test_helper.h:28
VERIFY
#define VERIFY(a)
Definition: main.h:380
EIGEN_TEST_MAX_SIZE
#define EIGEN_TEST_MAX_SIZE
Definition: boostmultiprec.cpp:16
EIGEN_DECLARE_TEST
EIGEN_DECLARE_TEST(linearstructure)
Definition: linearstructure.cpp:126
rows
int rows
Definition: Tutorial_commainit_02.cpp:1
real_complex
void real_complex(DenseIndex rows=MatrixType::RowsAtCompileTime, DenseIndex cols=MatrixType::ColsAtCompileTime)
Definition: linearstructure.cpp:76
CALL_SUBTEST_2
#define CALL_SUBTEST_2(FUNC)
Definition: split_test_helper.h:10
Eigen::Matrix
The matrix class, also used for vectors and row-vectors.
Definition: 3rdparty/Eigen/Eigen/src/Core/Matrix.h:178
abs
#define abs(x)
Definition: datatypes.h:17
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
linearstructure_overflow
void linearstructure_overflow()
Definition: linearstructure.cpp:114

gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:34:33