Functions
triangular.cpp File Reference
#include "main.h"

Functions

template<typename MatrixType >
void triangular_deprecated (const MatrixType &m)
 
template<typename MatrixType >
void triangular_square (const MatrixType &m)
 
template<typename MatrixType >
void triangular_rect (const MatrixType &m)
 
void bug_159 ()
 
 EIGEN_DECLARE_TEST (triangular)
 

Function Documentation

◆ bug_159()

void bug_159 ( )
261  {
262  Matrix3d m = Matrix3d::Random().triangularView<Lower>();
263  EIGEN_UNUSED_VARIABLE(m)
264 }
EIGEN_UNUSED_VARIABLE
#define EIGEN_UNUSED_VARIABLE(var)
Definition: Macros.h:966
Eigen::Lower
@ Lower
Definition: Constants.h:211
m
int * m
Definition: level2_cplx_impl.h:294

References EIGEN_UNUSED_VARIABLE, Eigen::Lower, and m.

Referenced by EIGEN_DECLARE_TEST().

◆ EIGEN_DECLARE_TEST()

EIGEN_DECLARE_TEST ( triangular  )
266  {
267  int maxsize = (std::min)(EIGEN_TEST_MAX_SIZE, 20);
268  for (int i = 0; i < g_repeat; i++) {
269  int r = internal::random<int>(2, maxsize);
270  TEST_SET_BUT_UNUSED_VARIABLE(r)
271  int c = internal::random<int>(2, maxsize);
272  TEST_SET_BUT_UNUSED_VARIABLE(c)
273 
274  CALL_SUBTEST_1(triangular_square(Matrix<float, 1, 1>()));
275  CALL_SUBTEST_2(triangular_square(Matrix<float, 2, 2>()));
276  CALL_SUBTEST_3(triangular_square(Matrix3d()));
277  CALL_SUBTEST_4(triangular_square(Matrix<std::complex<float>, 8, 8>()));
278  CALL_SUBTEST_5(triangular_square(MatrixXcd(r, r)));
279  CALL_SUBTEST_6(triangular_square(Matrix<float, Dynamic, Dynamic, RowMajor>(r, r)));
280 
281  CALL_SUBTEST_7(triangular_rect(Matrix<float, 4, 5>()));
282  CALL_SUBTEST_8(triangular_rect(Matrix<double, 6, 2>()));
283  CALL_SUBTEST_9(triangular_rect(MatrixXcf(r, c)));
284  CALL_SUBTEST_5(triangular_rect(MatrixXcd(r, c)));
285  CALL_SUBTEST_6(triangular_rect(Matrix<float, Dynamic, Dynamic, RowMajor>(r, c)));
286 
287  CALL_SUBTEST_100(triangular_deprecated(Matrix<float, 5, 7>()));
288  CALL_SUBTEST_100(triangular_deprecated(MatrixXd(r, c)));
289  }
290 
291  CALL_SUBTEST_1(bug_159());
292 }
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
maxsize
const size_t maxsize
Definition: benchmark-blocking-sizes.cpp:49
EIGEN_TEST_MAX_SIZE
#define EIGEN_TEST_MAX_SIZE
Definition: boostmultiprec.cpp:16
Eigen::Matrix
The matrix class, also used for vectors and row-vectors.
Definition: Eigen/Eigen/src/Core/Matrix.h:186
min
#define min(a, b)
Definition: datatypes.h:22
TEST_SET_BUT_UNUSED_VARIABLE
#define TEST_SET_BUT_UNUSED_VARIABLE(X)
Definition: main.h:139
Eigen::g_repeat
static int g_repeat
Definition: main.h:191
UniformPSDSelfTest.r
r
Definition: UniformPSDSelfTest.py:20
calibrate.c
int c
Definition: calibrate.py:100
CALL_SUBTEST_6
#define CALL_SUBTEST_6(FUNC)
Definition: split_test_helper.h:34
CALL_SUBTEST_3
#define CALL_SUBTEST_3(FUNC)
Definition: split_test_helper.h:16
CALL_SUBTEST_1
#define CALL_SUBTEST_1(FUNC)
Definition: split_test_helper.h:4
CALL_SUBTEST_8
#define CALL_SUBTEST_8(FUNC)
Definition: split_test_helper.h:46
CALL_SUBTEST_100
#define CALL_SUBTEST_100(FUNC)
Definition: split_test_helper.h:598
CALL_SUBTEST_5
#define CALL_SUBTEST_5(FUNC)
Definition: split_test_helper.h:28
CALL_SUBTEST_2
#define CALL_SUBTEST_2(FUNC)
Definition: split_test_helper.h:10
CALL_SUBTEST_7
#define CALL_SUBTEST_7(FUNC)
Definition: split_test_helper.h:40
CALL_SUBTEST_4
#define CALL_SUBTEST_4(FUNC)
Definition: split_test_helper.h:22
CALL_SUBTEST_9
#define CALL_SUBTEST_9(FUNC)
Definition: split_test_helper.h:52
triangular_deprecated
void triangular_deprecated(const MatrixType &m)
Definition: triangular.cpp:17
triangular_rect
void triangular_rect(const MatrixType &m)
Definition: triangular.cpp:187
triangular_square
void triangular_square(const MatrixType &m)
Definition: triangular.cpp:40
bug_159
void bug_159()
Definition: triangular.cpp:261

References bug_159(), calibrate::c, CALL_SUBTEST_1, CALL_SUBTEST_100, CALL_SUBTEST_2, CALL_SUBTEST_3, CALL_SUBTEST_4, CALL_SUBTEST_5, CALL_SUBTEST_6, CALL_SUBTEST_7, CALL_SUBTEST_8, CALL_SUBTEST_9, EIGEN_TEST_MAX_SIZE, Eigen::g_repeat, i, maxsize, min, UniformPSDSelfTest::r, TEST_SET_BUT_UNUSED_VARIABLE, triangular_deprecated(), triangular_rect(), and triangular_square().

◆ triangular_deprecated()

template<typename MatrixType >
void triangular_deprecated ( const MatrixType &  m)
17  {
18  Index rows = m.rows();
19  Index cols = m.cols();
20  MatrixType m1, m2, m3, m4;
21  m1.setRandom(rows, cols);
22  m2.setRandom(rows, cols);
23  m3 = m1;
24  m4 = m2;
25  // deprecated method:
26  m1.template triangularView<Eigen::Upper>().swap(m2);
27  // use this method instead:
28  m3.template triangularView<Eigen::Upper>().swap(m4.template triangularView<Eigen::Upper>());
29  VERIFY_IS_APPROX(m1, m3);
30  VERIFY_IS_APPROX(m2, m4);
31  // deprecated method:
32  m1.template triangularView<Eigen::Lower>().swap(m4);
33  // use this method instead:
34  m3.template triangularView<Eigen::Lower>().swap(m2.template triangularView<Eigen::Lower>());
35  VERIFY_IS_APPROX(m1, m3);
36  VERIFY_IS_APPROX(m2, m4);
37 }
m1
Matrix3d m1
Definition: IOFormat.cpp:2
m2
MatrixType m2(n_dims)
rows
int rows
Definition: Tutorial_commainit_02.cpp:1
cols
int cols
Definition: Tutorial_commainit_02.cpp:1
MatrixType
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
VERIFY_IS_APPROX
#define VERIFY_IS_APPROX(a, b)
Definition: integer_types.cpp:13
Eigen::Index
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:83

References cols, m, m1, m2(), rows, and VERIFY_IS_APPROX.

Referenced by EIGEN_DECLARE_TEST().

◆ triangular_rect()

template<typename MatrixType >
void triangular_rect ( const MatrixType &  m)
187  {
188  typedef typename MatrixType::Scalar Scalar;
189  typedef typename NumTraits<Scalar>::Real RealScalar;
190  enum { Rows = MatrixType::RowsAtCompileTime, Cols = MatrixType::ColsAtCompileTime };
191 
192  Index rows = m.rows();
193  Index cols = m.cols();
194 
195  MatrixType m1 = MatrixType::Random(rows, cols), m2 = MatrixType::Random(rows, cols), m3(rows, cols), m4(rows, cols),
196  r1(rows, cols), r2(rows, cols);
197 
198  MatrixType m1up = m1.template triangularView<Upper>();
199  MatrixType m2up = m2.template triangularView<Upper>();
200 
201  if (rows > 1 && cols > 1) {
202  VERIFY(m1up.isUpperTriangular());
203  VERIFY(m2up.transpose().isLowerTriangular());
204  VERIFY(!m2.isLowerTriangular());
205  }
206 
207  // test overloaded operator+=
208  r1.setZero();
209  r2.setZero();
210  r1.template triangularView<Upper>() += m1;
211  r2 += m1up;
212  VERIFY_IS_APPROX(r1, r2);
213 
214  // test overloaded operator=
215  m1.setZero();
216  m1.template triangularView<Upper>() = 3 * m2;
217  m3 = 3 * m2;
218  VERIFY_IS_APPROX(m3.template triangularView<Upper>().toDenseMatrix(), m1);
219 
220  m1.setZero();
221  m1.template triangularView<Lower>() = 3 * m2;
222  VERIFY_IS_APPROX(m3.template triangularView<Lower>().toDenseMatrix(), m1);
223 
224  m1.setZero();
225  m1.template triangularView<StrictlyUpper>() = 3 * m2;
226  VERIFY_IS_APPROX(m3.template triangularView<StrictlyUpper>().toDenseMatrix(), m1);
227 
228  m1.setZero();
229  m1.template triangularView<StrictlyLower>() = 3 * m2;
230  VERIFY_IS_APPROX(m3.template triangularView<StrictlyLower>().toDenseMatrix(), m1);
231  m1.setRandom();
232  m2 = m1.template triangularView<Upper>();
233  VERIFY(m2.isUpperTriangular());
234  VERIFY(!m2.isLowerTriangular());
235  m2 = m1.template triangularView<StrictlyUpper>();
236  VERIFY(m2.isUpperTriangular());
237  VERIFY(m2.diagonal().isMuchSmallerThan(RealScalar(1)));
238  m2 = m1.template triangularView<UnitUpper>();
239  VERIFY(m2.isUpperTriangular());
240  m2.diagonal().array() -= Scalar(1);
241  VERIFY(m2.diagonal().isMuchSmallerThan(RealScalar(1)));
242  m2 = m1.template triangularView<Lower>();
243  VERIFY(m2.isLowerTriangular());
244  VERIFY(!m2.isUpperTriangular());
245  m2 = m1.template triangularView<StrictlyLower>();
246  VERIFY(m2.isLowerTriangular());
247  VERIFY(m2.diagonal().isMuchSmallerThan(RealScalar(1)));
248  m2 = m1.template triangularView<UnitLower>();
249  VERIFY(m2.isLowerTriangular());
250  m2.diagonal().array() -= Scalar(1);
251  VERIFY(m2.diagonal().isMuchSmallerThan(RealScalar(1)));
252  // test swap
253  m1.setOnes();
254  m2.setZero();
255  m2.template triangularView<Upper>().swap(m1.template triangularView<Eigen::Upper>());
256  m3.setZero();
257  m3.template triangularView<Upper>().setOnes();
258  VERIFY_IS_APPROX(m2, m3);
259 }
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:45
RealScalar
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
VERIFY
#define VERIFY(a)
Definition: main.h:362
Eigen::NumTraits
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:217

References cols, m, m1, m2(), rows, VERIFY, and VERIFY_IS_APPROX.

Referenced by EIGEN_DECLARE_TEST().

◆ triangular_square()

template<typename MatrixType >
void triangular_square ( const MatrixType &  m)
40  {
41  typedef typename MatrixType::Scalar Scalar;
42  typedef typename NumTraits<Scalar>::Real RealScalar;
43  typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType;
44 
45  RealScalar largerEps = 10 * test_precision<RealScalar>();
46 
47  Index rows = m.rows();
48  Index cols = m.cols();
49 
50  MatrixType m1 = MatrixType::Random(rows, cols), m2 = MatrixType::Random(rows, cols), m3(rows, cols), m4(rows, cols),
51  r1(rows, cols), r2(rows, cols);
52  VectorType v2 = VectorType::Random(rows);
53  VectorType v3 = VectorType::Zero(rows);
54 
55  MatrixType m1up = m1.template triangularView<Upper>();
56  MatrixType m2up = m2.template triangularView<Upper>();
57 
58  if (rows * cols > 1) {
59  VERIFY(m1up.isUpperTriangular());
60  VERIFY(m2up.transpose().isLowerTriangular());
61  VERIFY(!m2.isLowerTriangular());
62  }
63 
64  // VERIFY_IS_APPROX(m1up.transpose() * m2, m1.upper().transpose().lower() * m2);
65 
66  // test overloaded operator+=
67  r1.setZero();
68  r2.setZero();
69  r1.template triangularView<Upper>() += m1;
70  r2 += m1up;
71  VERIFY_IS_APPROX(r1, r2);
72 
73  // test overloaded operator=
74  m1.setZero();
75  m1.template triangularView<Upper>() = m2.transpose() + m2;
76  m3 = m2.transpose() + m2;
77  VERIFY_IS_APPROX(m3.template triangularView<Lower>().transpose().toDenseMatrix(), m1);
78 
79  // test overloaded operator=
80  m1.setZero();
81  m1.template triangularView<Lower>() = m2.transpose() + m2;
82  VERIFY_IS_APPROX(m3.template triangularView<Lower>().toDenseMatrix(), m1);
83 
84  VERIFY_IS_APPROX(m3.template triangularView<Lower>().conjugate().toDenseMatrix(),
85  m3.conjugate().template triangularView<Lower>().toDenseMatrix());
86 
87  m1 = MatrixType::Random(rows, cols);
88  for (int i = 0; i < rows; ++i)
89  while (numext::abs2(m1(i, i)) < RealScalar(1e-1)) m1(i, i) = internal::random<Scalar>();
90 
91  Transpose<MatrixType> trm4(m4);
92  // test back and forward substitution with a vector as the rhs
93  m3 = m1.template triangularView<Upper>();
94  v3 = m3.adjoint() * (m1.adjoint().template triangularView<Lower>().solve(v2));
95  VERIFY(v2.isApprox(v3, largerEps));
96  m3 = m1.template triangularView<Lower>();
97  v3 = m3.transpose() * (m1.transpose().template triangularView<Upper>().solve(v2));
98  VERIFY(v2.isApprox(v3, largerEps));
99  m3 = m1.template triangularView<Upper>();
100  v3 = m3 * (m1.template triangularView<Upper>().solve(v2));
101  VERIFY(v2.isApprox(v3, largerEps));
102  m3 = m1.template triangularView<Lower>();
103  v3 = m3.conjugate() * (m1.conjugate().template triangularView<Lower>().solve(v2));
104  VERIFY(v2.isApprox(v3, largerEps));
105 
106  // test back and forward substitution with a matrix as the rhs
107  m3 = m1.template triangularView<Upper>();
108  m4 = m3.adjoint() * (m1.adjoint().template triangularView<Lower>().solve(m2));
109  VERIFY(m2.isApprox(m4, largerEps));
110  m3 = m1.template triangularView<Lower>();
111  m4 = m3.transpose() * (m1.transpose().template triangularView<Upper>().solve(m2));
112  VERIFY(m2.isApprox(m4, largerEps));
113  m3 = m1.template triangularView<Upper>();
114  m4 = m3 * (m1.template triangularView<Upper>().solve(m2));
115  VERIFY(m2.isApprox(m4, largerEps));
116  m3 = m1.template triangularView<Lower>();
117  m4 = m3.conjugate() * (m1.conjugate().template triangularView<Lower>().solve(m2));
118  VERIFY(m2.isApprox(m4, largerEps));
119 
120  // check M * inv(L) using in place API
121  m4 = m3;
122  m1.transpose().template triangularView<Eigen::Upper>().solveInPlace(trm4);
123  VERIFY_IS_APPROX(m4 * m1.template triangularView<Eigen::Lower>(), m3);
124 
125  // check M * inv(U) using in place API
126  m3 = m1.template triangularView<Upper>();
127  m4 = m3;
128  m3.transpose().template triangularView<Eigen::Lower>().solveInPlace(trm4);
129  VERIFY_IS_APPROX(m4 * m1.template triangularView<Eigen::Upper>(), m3);
130 
131  // check solve with unit diagonal
132  m3 = m1.template triangularView<UnitUpper>();
133  VERIFY(m2.isApprox(m3 * (m1.template triangularView<UnitUpper>().solve(m2)), largerEps));
134 
135  // VERIFY(( m1.template triangularView<Upper>()
136  // * m2.template triangularView<Upper>()).isUpperTriangular());
137 
138  // test swap
139  m1.setOnes();
140  m2.setZero();
141  m2.template triangularView<Upper>().swap(m1.template triangularView<Eigen::Upper>());
142  m3.setZero();
143  m3.template triangularView<Upper>().setOnes();
144  VERIFY_IS_APPROX(m2, m3);
145 
146  m1.setRandom();
147  m3 = m1.template triangularView<Upper>();
148  Matrix<Scalar, MatrixType::ColsAtCompileTime, Dynamic> m5(cols, internal::random<int>(1, 20));
149  m5.setRandom();
150  Matrix<Scalar, Dynamic, MatrixType::RowsAtCompileTime> m6(internal::random<int>(1, 20), rows);
151  m6.setRandom();
152  VERIFY_IS_APPROX(m1.template triangularView<Upper>() * m5, m3 * m5);
153  VERIFY_IS_APPROX(m6 * m1.template triangularView<Upper>(), m6 * m3);
154 
155  m1up = m1.template triangularView<Upper>();
156  VERIFY_IS_APPROX(m1.template selfadjointView<Upper>().template triangularView<Upper>().toDenseMatrix(), m1up);
157  VERIFY_IS_APPROX(m1up.template selfadjointView<Upper>().template triangularView<Upper>().toDenseMatrix(), m1up);
158  VERIFY_IS_APPROX(m1.template selfadjointView<Upper>().template triangularView<Lower>().toDenseMatrix(),
159  m1up.adjoint());
160  VERIFY_IS_APPROX(m1up.template selfadjointView<Upper>().template triangularView<Lower>().toDenseMatrix(),
161  m1up.adjoint());
162 
163  VERIFY_IS_APPROX(m1.template selfadjointView<Upper>().diagonal(), m1.diagonal());
164 
165  m3.setRandom();
166  const MatrixType& m3c(m3);
167  VERIFY(is_same_type(m3c.template triangularView<Lower>(),
168  m3.template triangularView<Lower>().template conjugateIf<false>()));
169  VERIFY(is_same_type(m3c.template triangularView<Lower>().conjugate(),
170  m3.template triangularView<Lower>().template conjugateIf<true>()));
171  VERIFY_IS_APPROX(m3.template triangularView<Lower>().template conjugateIf<true>().toDenseMatrix(),
172  m3.conjugate().template triangularView<Lower>().toDenseMatrix());
173  VERIFY_IS_APPROX(m3.template triangularView<Lower>().template conjugateIf<false>().toDenseMatrix(),
174  m3.template triangularView<Lower>().toDenseMatrix());
175 
176  VERIFY(is_same_type(m3c.template selfadjointView<Lower>(),
177  m3.template selfadjointView<Lower>().template conjugateIf<false>()));
178  VERIFY(is_same_type(m3c.template selfadjointView<Lower>().conjugate(),
179  m3.template selfadjointView<Lower>().template conjugateIf<true>()));
180  VERIFY_IS_APPROX(m3.template selfadjointView<Lower>().template conjugateIf<true>().toDenseMatrix(),
181  m3.conjugate().template selfadjointView<Lower>().toDenseMatrix());
182  VERIFY_IS_APPROX(m3.template selfadjointView<Lower>().template conjugateIf<false>().toDenseMatrix(),
183  m3.template selfadjointView<Lower>().toDenseMatrix());
184 }
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
v2
Map< RowVectorXf > v2(M2.data(), M2.size())
Eigen::Transpose
Expression of the transpose of a matrix.
Definition: Transpose.h:56
Eigen::numext::abs2
EIGEN_DEVICE_FUNC bool abs2(bool x)
Definition: MathFunctions.h:1102
Eigen::is_same_type
std::enable_if_t< internal::is_same< T1, T2 >::value, bool > is_same_type(const T1 &, const T2 &)
Definition: main.h:407
oomph::PseudoSolidHelper::Zero
double Zero
Definition: pseudosolid_node_update_elements.cc:35
VectorType
Definition: fft_test_shared.h:66

References Eigen::numext::abs2(), cols, e(), i, Eigen::is_same_type(), m, m1, m2(), rows, Eigen::PlainObjectBase< Derived >::setRandom(), v2(), VERIFY, VERIFY_IS_APPROX, and oomph::PseudoSolidHelper::Zero.

Referenced by EIGEN_DECLARE_TEST().