1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2011 Benoit Jacob <jacob.benoit.1@gmail.com>
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#define EIGEN_NO_STATIC_ASSERT
11
12#include "main.h"
13
14template<typename ArrayType> void vectorwiseop_array(const ArrayType& m)
15{
16  typedef typename ArrayType::Index Index;
17  typedef typename ArrayType::Scalar Scalar;
18  typedef Array<Scalar, ArrayType::RowsAtCompileTime, 1> ColVectorType;
19  typedef Array<Scalar, 1, ArrayType::ColsAtCompileTime> RowVectorType;
20
21  Index rows = m.rows();
22  Index cols = m.cols();
23  Index r = internal::random<Index>(0, rows-1),
24        c = internal::random<Index>(0, cols-1);
25
26  ArrayType m1 = ArrayType::Random(rows, cols),
27            m2(rows, cols),
28            m3(rows, cols);
29
30  ColVectorType colvec = ColVectorType::Random(rows);
31  RowVectorType rowvec = RowVectorType::Random(cols);
32
33  // test addition
34
35  m2 = m1;
36  m2.colwise() += colvec;
37  VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
38  VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
39
40  VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
41  VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
42
43  m2 = m1;
44  m2.rowwise() += rowvec;
45  VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
46  VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
47
48  VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
49  VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
50
51  // test substraction
52
53  m2 = m1;
54  m2.colwise() -= colvec;
55  VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
56  VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
57
58  VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
59  VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
60
61  m2 = m1;
62  m2.rowwise() -= rowvec;
63  VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
64  VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
65
66  VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
67  VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
68
69  // test multiplication
70
71  m2 = m1;
72  m2.colwise() *= colvec;
73  VERIFY_IS_APPROX(m2, m1.colwise() * colvec);
74  VERIFY_IS_APPROX(m2.col(c), m1.col(c) * colvec);
75
76  VERIFY_RAISES_ASSERT(m2.colwise() *= colvec.transpose());
77  VERIFY_RAISES_ASSERT(m1.colwise() * colvec.transpose());
78
79  m2 = m1;
80  m2.rowwise() *= rowvec;
81  VERIFY_IS_APPROX(m2, m1.rowwise() * rowvec);
82  VERIFY_IS_APPROX(m2.row(r), m1.row(r) * rowvec);
83
84  VERIFY_RAISES_ASSERT(m2.rowwise() *= rowvec.transpose());
85  VERIFY_RAISES_ASSERT(m1.rowwise() * rowvec.transpose());
86
87  // test quotient
88
89  m2 = m1;
90  m2.colwise() /= colvec;
91  VERIFY_IS_APPROX(m2, m1.colwise() / colvec);
92  VERIFY_IS_APPROX(m2.col(c), m1.col(c) / colvec);
93
94  VERIFY_RAISES_ASSERT(m2.colwise() /= colvec.transpose());
95  VERIFY_RAISES_ASSERT(m1.colwise() / colvec.transpose());
96
97  m2 = m1;
98  m2.rowwise() /= rowvec;
99  VERIFY_IS_APPROX(m2, m1.rowwise() / rowvec);
100  VERIFY_IS_APPROX(m2.row(r), m1.row(r) / rowvec);
101
102  VERIFY_RAISES_ASSERT(m2.rowwise() /= rowvec.transpose());
103  VERIFY_RAISES_ASSERT(m1.rowwise() / rowvec.transpose());
104
105  m2 = m1;
106  // yes, there might be an aliasing issue there but ".rowwise() /="
107  // is suppposed to evaluate " m2.colwise().sum()" into to temporary to avoid
108  // evaluating the reducions multiple times
109  if(ArrayType::RowsAtCompileTime>2 || ArrayType::RowsAtCompileTime==Dynamic)
110  {
111    m2.rowwise() /= m2.colwise().sum();
112    VERIFY_IS_APPROX(m2, m1.rowwise() / m1.colwise().sum());
113  }
114}
115
116template<typename MatrixType> void vectorwiseop_matrix(const MatrixType& m)
117{
118  typedef typename MatrixType::Index Index;
119  typedef typename MatrixType::Scalar Scalar;
120  typedef typename NumTraits<Scalar>::Real RealScalar;
121  typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType;
122  typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType;
123  typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, 1> RealColVectorType;
124  typedef Matrix<RealScalar, 1, MatrixType::ColsAtCompileTime> RealRowVectorType;
125
126  Index rows = m.rows();
127  Index cols = m.cols();
128  Index r = internal::random<Index>(0, rows-1),
129        c = internal::random<Index>(0, cols-1);
130
131  MatrixType m1 = MatrixType::Random(rows, cols),
132            m2(rows, cols),
133            m3(rows, cols);
134
135  ColVectorType colvec = ColVectorType::Random(rows);
136  RowVectorType rowvec = RowVectorType::Random(cols);
137  RealColVectorType rcres;
138  RealRowVectorType rrres;
139
140  // test addition
141
142  m2 = m1;
143  m2.colwise() += colvec;
144  VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
145  VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
146
147  VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
148  VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
149
150  m2 = m1;
151  m2.rowwise() += rowvec;
152  VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
153  VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
154
155  VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
156  VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
157
158  // test substraction
159
160  m2 = m1;
161  m2.colwise() -= colvec;
162  VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
163  VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
164
165  VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
166  VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
167
168  m2 = m1;
169  m2.rowwise() -= rowvec;
170  VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
171  VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
172
173  VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
174  VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
175
176  // test norm
177  rrres = m1.colwise().norm();
178  VERIFY_IS_APPROX(rrres(c), m1.col(c).norm());
179  rcres = m1.rowwise().norm();
180  VERIFY_IS_APPROX(rcres(r), m1.row(r).norm());
181
182  // test normalized
183  m2 = m1.colwise().normalized();
184  VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
185  m2 = m1.rowwise().normalized();
186  VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
187
188  // test normalize
189  m2 = m1;
190  m2.colwise().normalize();
191  VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
192  m2 = m1;
193  m2.rowwise().normalize();
194  VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
195}
196
197void test_vectorwiseop()
198{
199  CALL_SUBTEST_1(vectorwiseop_array(Array22cd()));
200  CALL_SUBTEST_2(vectorwiseop_array(Array<double, 3, 2>()));
201  CALL_SUBTEST_3(vectorwiseop_array(ArrayXXf(3, 4)));
202  CALL_SUBTEST_4(vectorwiseop_matrix(Matrix4cf()));
203  CALL_SUBTEST_5(vectorwiseop_matrix(Matrix<float,4,5>()));
204  CALL_SUBTEST_6(vectorwiseop_matrix(MatrixXd(7,2)));
205}
206