product_trmv.cpp
Go to the documentation of this file.
1 // This file is triangularView of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2009 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 #include "main.h"
11 
12 template<typename MatrixType> void trmv(const MatrixType& m)
13 {
14  typedef typename MatrixType::Scalar Scalar;
15  typedef typename NumTraits<Scalar>::Real RealScalar;
17 
18  RealScalar largerEps = 10*test_precision<RealScalar>();
19 
20  Index rows = m.rows();
21  Index cols = m.cols();
22 
23  MatrixType m1 = MatrixType::Random(rows, cols),
24  m3(rows, cols);
25  VectorType v1 = VectorType::Random(rows);
26 
27  Scalar s1 = internal::random<Scalar>();
28 
29  m1 = MatrixType::Random(rows, cols);
30 
31  // check with a column-major matrix
32  m3 = m1.template triangularView<Eigen::Lower>();
33  VERIFY((m3 * v1).isApprox(m1.template triangularView<Eigen::Lower>() * v1, largerEps));
34  m3 = m1.template triangularView<Eigen::Upper>();
35  VERIFY((m3 * v1).isApprox(m1.template triangularView<Eigen::Upper>() * v1, largerEps));
36  m3 = m1.template triangularView<Eigen::UnitLower>();
37  VERIFY((m3 * v1).isApprox(m1.template triangularView<Eigen::UnitLower>() * v1, largerEps));
38  m3 = m1.template triangularView<Eigen::UnitUpper>();
39  VERIFY((m3 * v1).isApprox(m1.template triangularView<Eigen::UnitUpper>() * v1, largerEps));
40 
41  // check conjugated and scalar multiple expressions (col-major)
42  m3 = m1.template triangularView<Eigen::Lower>();
43  VERIFY(((s1*m3).conjugate() * v1).isApprox((s1*m1).conjugate().template triangularView<Eigen::Lower>() * v1, largerEps));
44  m3 = m1.template triangularView<Eigen::Upper>();
45  VERIFY((m3.conjugate() * v1.conjugate()).isApprox(m1.conjugate().template triangularView<Eigen::Upper>() * v1.conjugate(), largerEps));
46 
47  // check with a row-major matrix
48  m3 = m1.template triangularView<Eigen::Upper>();
49  VERIFY((m3.transpose() * v1).isApprox(m1.transpose().template triangularView<Eigen::Lower>() * v1, largerEps));
50  m3 = m1.template triangularView<Eigen::Lower>();
51  VERIFY((m3.transpose() * v1).isApprox(m1.transpose().template triangularView<Eigen::Upper>() * v1, largerEps));
52  m3 = m1.template triangularView<Eigen::UnitUpper>();
53  VERIFY((m3.transpose() * v1).isApprox(m1.transpose().template triangularView<Eigen::UnitLower>() * v1, largerEps));
54  m3 = m1.template triangularView<Eigen::UnitLower>();
55  VERIFY((m3.transpose() * v1).isApprox(m1.transpose().template triangularView<Eigen::UnitUpper>() * v1, largerEps));
56 
57  // check conjugated and scalar multiple expressions (row-major)
58  m3 = m1.template triangularView<Eigen::Upper>();
59  VERIFY((m3.adjoint() * v1).isApprox(m1.adjoint().template triangularView<Eigen::Lower>() * v1, largerEps));
60  m3 = m1.template triangularView<Eigen::Lower>();
61  VERIFY((m3.adjoint() * (s1*v1.conjugate())).isApprox(m1.adjoint().template triangularView<Eigen::Upper>() * (s1*v1.conjugate()), largerEps));
62  m3 = m1.template triangularView<Eigen::UnitUpper>();
63 
64  // check transposed cases:
65  m3 = m1.template triangularView<Eigen::Lower>();
66  VERIFY((v1.transpose() * m3).isApprox(v1.transpose() * m1.template triangularView<Eigen::Lower>(), largerEps));
67  VERIFY((v1.adjoint() * m3).isApprox(v1.adjoint() * m1.template triangularView<Eigen::Lower>(), largerEps));
68  VERIFY((v1.adjoint() * m3.adjoint()).isApprox(v1.adjoint() * m1.template triangularView<Eigen::Lower>().adjoint(), largerEps));
69 
70  // TODO check with sub-matrices
71 }
72 
73 EIGEN_DECLARE_TEST(product_trmv)
74 {
75  int s = 0;
76  for(int i = 0; i < g_repeat ; i++) {
79  CALL_SUBTEST_3( trmv(Matrix3d()) );
80 
81  s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2);
82  CALL_SUBTEST_4( trmv(MatrixXcf(s,s)) );
83  CALL_SUBTEST_5( trmv(MatrixXcd(s,s)) );
85 
86  s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE);
89  }
90 }
Matrix3f m
SCALAR Scalar
Definition: bench_gemm.cpp:46
#define CALL_SUBTEST_6(FUNC)
#define CALL_SUBTEST_4(FUNC)
void adjoint(const MatrixType &m)
Definition: adjoint.cpp:67
Vector v1
#define CALL_SUBTEST_3(FUNC)
MatrixXf MatrixType
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
EIGEN_DECLARE_TEST(product_trmv)
void trmv(const MatrixType &m)
#define CALL_SUBTEST_1(FUNC)
Matrix3d m1
Definition: IOFormat.cpp:2
static int g_repeat
Definition: main.h:169
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
RealScalar s
EIGEN_DEVICE_FUNC ConjugateReturnType conjugate() const
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
#define TEST_SET_BUT_UNUSED_VARIABLE(X)
Definition: main.h:121
#define CALL_SUBTEST_5(FUNC)
#define VERIFY(a)
Definition: main.h:380
#define EIGEN_TEST_MAX_SIZE
#define CALL_SUBTEST_2(FUNC)
The matrix class, also used for vectors and row-vectors.
EIGEN_DEVICE_FUNC bool isApprox(const Scalar &x, const Scalar &y, const typename NumTraits< Scalar >::Real &precision=NumTraits< Scalar >::dummy_precision())


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:35:19