Functions
sparse_block.cpp File Reference
#include "sparse.h"
#include "AnnoyingScalar.h"

Functions

template<typename T >
std::enable_if_t<(T::Flags &RowMajorBit)==RowMajorBit, typename T::RowXpr > innervec (T &A, Index i)
 
template<typename T >
std::enable_if_t<(T::Flags &RowMajorBit)==0, typename T::ColXpr > innervec (T &A, Index i)
 
template<typename SparseMatrixType >
void sparse_block (const SparseMatrixType &ref)
 
 EIGEN_DECLARE_TEST (sparse_block)
 

Function Documentation

◆ EIGEN_DECLARE_TEST()

EIGEN_DECLARE_TEST ( sparse_block  )
285  {
286  for (int i = 0; i < g_repeat; i++) {
287  int r = Eigen::internal::random<int>(1, 200), c = Eigen::internal::random<int>(1, 200);
288  if (Eigen::internal::random<int>(0, 4) == 0) {
289  r = c; // check square matrices in 25% of tries
290  }
291  EIGEN_UNUSED_VARIABLE(r + c);
292  CALL_SUBTEST_1((sparse_block(SparseMatrix<double>(1, 1))));
293  CALL_SUBTEST_1((sparse_block(SparseMatrix<double>(8, 8))));
294  CALL_SUBTEST_1((sparse_block(SparseMatrix<double>(r, c))));
295  CALL_SUBTEST_2((sparse_block(SparseMatrix<std::complex<double>, ColMajor>(r, c))));
296  CALL_SUBTEST_2((sparse_block(SparseMatrix<std::complex<double>, RowMajor>(r, c))));
297 
298  CALL_SUBTEST_3((sparse_block(SparseMatrix<double, ColMajor, long int>(r, c))));
299  CALL_SUBTEST_3((sparse_block(SparseMatrix<double, RowMajor, long int>(r, c))));
300 
301  r = Eigen::internal::random<int>(1, 100);
302  c = Eigen::internal::random<int>(1, 100);
303  if (Eigen::internal::random<int>(0, 4) == 0) {
304  r = c; // check square matrices in 25% of tries
305  }
306 
307  CALL_SUBTEST_4((sparse_block(SparseMatrix<double, ColMajor, short int>(short(r), short(c)))));
308  CALL_SUBTEST_4((sparse_block(SparseMatrix<double, RowMajor, short int>(short(r), short(c)))));
309 #ifndef EIGEN_TEST_ANNOYING_SCALAR_DONT_THROW
310  AnnoyingScalar::dont_throw = true;
311 #endif
312  CALL_SUBTEST_5((sparse_block(SparseMatrix<AnnoyingScalar>(r, c))));
313  }
314 }
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
EIGEN_UNUSED_VARIABLE
#define EIGEN_UNUSED_VARIABLE(var)
Definition: Macros.h:966
AnnoyingScalar::dont_throw
static bool dont_throw
Definition: AnnoyingScalar.h:127
Eigen::ColMajor
@ ColMajor
Definition: Constants.h:318
Eigen::RowMajor
@ RowMajor
Definition: Constants.h:320
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
sparse_block
void sparse_block(const SparseMatrixType &ref)
Definition: sparse_block.cpp:24
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_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_4
#define CALL_SUBTEST_4(FUNC)
Definition: split_test_helper.h:22

References calibrate::c, CALL_SUBTEST_1, CALL_SUBTEST_2, CALL_SUBTEST_3, CALL_SUBTEST_4, CALL_SUBTEST_5, Eigen::ColMajor, AnnoyingScalar::dont_throw, EIGEN_UNUSED_VARIABLE, Eigen::g_repeat, i, UniformPSDSelfTest::r, Eigen::RowMajor, and sparse_block().

◆ innervec() [1/2]

template<typename T >
std::enable_if_t<(T::Flags & RowMajorBit) == RowMajorBit, typename T::RowXpr> innervec ( T &  A,
Index  i 
)
14  {
15  return A.row(i);
16 }
Eigen::Matrix
The matrix class, also used for vectors and row-vectors.
Definition: Eigen/Eigen/src/Core/Matrix.h:186

References i.

Referenced by sparse_block().

◆ innervec() [2/2]

template<typename T >
std::enable_if_t<(T::Flags & RowMajorBit) == 0, typename T::ColXpr> innervec ( T &  A,
Index  i 
)
19  {
20  return A.col(i);
21 }

References i.

◆ sparse_block()

template<typename SparseMatrixType >
void sparse_block ( const SparseMatrixType &  ref)
24  {
25  const Index rows = ref.rows();
26  const Index cols = ref.cols();
27  const Index inner = ref.innerSize();
28  const Index outer = ref.outerSize();
29 
30  typedef typename SparseMatrixType::Scalar Scalar;
31  typedef typename SparseMatrixType::RealScalar RealScalar;
32  typedef typename SparseMatrixType::StorageIndex StorageIndex;
33 
34  double density = (std::max)(8. / (rows * cols), 0.01);
35  typedef Matrix<Scalar, Dynamic, Dynamic, SparseMatrixType::IsRowMajor ? RowMajor : ColMajor> DenseMatrix;
36  typedef Matrix<Scalar, Dynamic, 1> DenseVector;
37  typedef Matrix<Scalar, 1, Dynamic> RowDenseVector;
38  typedef SparseVector<Scalar> SparseVectorType;
39 
40  Scalar s1 = internal::random<Scalar>();
41  {
42  SparseMatrixType m(rows, cols);
43  DenseMatrix refMat = DenseMatrix::Zero(rows, cols);
44  initSparse<Scalar>(density, refMat, m);
45 
46  VERIFY_IS_APPROX(m, refMat);
47 
48  // test InnerIterators and Block expressions
49  for (int t = 0; t < 10; ++t) {
50  Index j = internal::random<Index>(0, cols - 2);
51  Index i = internal::random<Index>(0, rows - 2);
52  Index w = internal::random<Index>(1, cols - j);
53  Index h = internal::random<Index>(1, rows - i);
54 
55  VERIFY_IS_APPROX(m.block(i, j, h, w), refMat.block(i, j, h, w));
56  for (Index c = 0; c < w; c++) {
57  VERIFY_IS_APPROX(m.block(i, j, h, w).col(c), refMat.block(i, j, h, w).col(c));
58  for (Index r = 0; r < h; r++) {
59  VERIFY_IS_APPROX(m.block(i, j, h, w).col(c).coeff(r), refMat.block(i, j, h, w).col(c).coeff(r));
60  VERIFY_IS_APPROX(m.block(i, j, h, w).coeff(r, c), refMat.block(i, j, h, w).coeff(r, c));
61  }
62  }
63  for (Index r = 0; r < h; r++) {
64  VERIFY_IS_APPROX(m.block(i, j, h, w).row(r), refMat.block(i, j, h, w).row(r));
65  for (Index c = 0; c < w; c++) {
66  VERIFY_IS_APPROX(m.block(i, j, h, w).row(r).coeff(c), refMat.block(i, j, h, w).row(r).coeff(c));
67  VERIFY_IS_APPROX(m.block(i, j, h, w).coeff(r, c), refMat.block(i, j, h, w).coeff(r, c));
68  }
69  }
70 
71  VERIFY_IS_APPROX(m.middleCols(j, w), refMat.middleCols(j, w));
72  VERIFY_IS_APPROX(m.middleRows(i, h), refMat.middleRows(i, h));
73  for (Index r = 0; r < h; r++) {
74  VERIFY_IS_APPROX(m.middleCols(j, w).row(r), refMat.middleCols(j, w).row(r));
75  VERIFY_IS_APPROX(m.middleRows(i, h).row(r), refMat.middleRows(i, h).row(r));
76  for (Index c = 0; c < w; c++) {
77  VERIFY_IS_APPROX(m.col(c).coeff(r), refMat.col(c).coeff(r));
78  VERIFY_IS_APPROX(m.row(r).coeff(c), refMat.row(r).coeff(c));
79 
80  VERIFY_IS_APPROX(m.middleCols(j, w).coeff(r, c), refMat.middleCols(j, w).coeff(r, c));
81  VERIFY_IS_APPROX(m.middleRows(i, h).coeff(r, c), refMat.middleRows(i, h).coeff(r, c));
82  if (!numext::is_exactly_zero(m.middleCols(j, w).coeff(r, c))) {
83  VERIFY_IS_APPROX(m.middleCols(j, w).coeffRef(r, c), refMat.middleCols(j, w).coeff(r, c));
84  }
85  if (!numext::is_exactly_zero(m.middleRows(i, h).coeff(r, c))) {
86  VERIFY_IS_APPROX(m.middleRows(i, h).coeff(r, c), refMat.middleRows(i, h).coeff(r, c));
87  }
88  }
89  }
90  for (Index c = 0; c < w; c++) {
91  VERIFY_IS_APPROX(m.middleCols(j, w).col(c), refMat.middleCols(j, w).col(c));
92  VERIFY_IS_APPROX(m.middleRows(i, h).col(c), refMat.middleRows(i, h).col(c));
93  }
94  }
95 
96  for (Index c = 0; c < cols; c++) {
97  VERIFY_IS_APPROX(m.col(c) + m.col(c), (m + m).col(c));
98  VERIFY_IS_APPROX(m.col(c) + m.col(c), refMat.col(c) + refMat.col(c));
99  }
100 
101  for (Index r = 0; r < rows; r++) {
102  VERIFY_IS_APPROX(m.row(r) + m.row(r), (m + m).row(r));
103  VERIFY_IS_APPROX(m.row(r) + m.row(r), refMat.row(r) + refMat.row(r));
104  }
105  }
106 
107  // test innerVector()
108  {
109  DenseMatrix refMat2 = DenseMatrix::Zero(rows, cols);
110  SparseMatrixType m2(rows, cols);
111  initSparse<Scalar>(density, refMat2, m2);
112  Index j0 = internal::random<Index>(0, outer - 1);
113  Index j1 = internal::random<Index>(0, outer - 1);
114  Index r0 = internal::random<Index>(0, rows - 1);
115  Index c0 = internal::random<Index>(0, cols - 1);
116 
117  VERIFY_IS_APPROX(m2.innerVector(j0), innervec(refMat2, j0));
118  VERIFY_IS_APPROX(m2.innerVector(j0) + m2.innerVector(j1), innervec(refMat2, j0) + innervec(refMat2, j1));
119 
120  m2.innerVector(j0) *= Scalar(2);
121  innervec(refMat2, j0) *= Scalar(2);
122  VERIFY_IS_APPROX(m2, refMat2);
123 
124  m2.row(r0) *= Scalar(3);
125  refMat2.row(r0) *= Scalar(3);
126  VERIFY_IS_APPROX(m2, refMat2);
127 
128  m2.col(c0) *= Scalar(4);
129  refMat2.col(c0) *= Scalar(4);
130  VERIFY_IS_APPROX(m2, refMat2);
131 
132  m2.row(r0) /= Scalar(3);
133  refMat2.row(r0) /= Scalar(3);
134  VERIFY_IS_APPROX(m2, refMat2);
135 
136  m2.col(c0) /= Scalar(4);
137  refMat2.col(c0) /= Scalar(4);
138  VERIFY_IS_APPROX(m2, refMat2);
139 
140  SparseVectorType v1;
141  VERIFY_IS_APPROX(v1 = m2.col(c0) * 4, refMat2.col(c0) * 4);
142  VERIFY_IS_APPROX(v1 = m2.row(r0) * 4, refMat2.row(r0).transpose() * 4);
143 
144  SparseMatrixType m3(rows, cols);
145  m3.reserve(VectorXi::Constant(outer, int(inner / 2)));
146  for (Index j = 0; j < outer; ++j)
147  for (Index k = 0; k < (std::min)(j, inner); ++k)
148  m3.insertByOuterInner(j, k) = internal::convert_index<StorageIndex>(k + 1);
149  for (Index j = 0; j < (std::min)(outer, inner); ++j) {
150  VERIFY(j == numext::real(m3.innerVector(j).nonZeros()));
151  if (j > 0) VERIFY_IS_EQUAL(RealScalar(j), numext::real(m3.innerVector(j).lastCoeff()));
152  }
153  m3.makeCompressed();
154  for (Index j = 0; j < (std::min)(outer, inner); ++j) {
155  VERIFY(j == numext::real(m3.innerVector(j).nonZeros()));
156  if (j > 0) VERIFY_IS_EQUAL(RealScalar(j), numext::real(m3.innerVector(j).lastCoeff()));
157  }
158 
159  VERIFY(m3.innerVector(j0).nonZeros() == m3.transpose().innerVector(j0).nonZeros());
160 
161  // m2.innerVector(j0) = 2*m2.innerVector(j1);
162  // refMat2.col(j0) = 2*refMat2.col(j1);
163  // VERIFY_IS_APPROX(m2, refMat2);
164  }
165 
166  // test innerVectors()
167  {
168  DenseMatrix refMat2 = DenseMatrix::Zero(rows, cols);
169  SparseMatrixType m2(rows, cols);
170  initSparse<Scalar>(density, refMat2, m2);
171  if (internal::random<float>(0, 1) > 0.5f) m2.makeCompressed();
172  Index j0 = internal::random<Index>(0, outer - 2);
173  Index j1 = internal::random<Index>(0, outer - 2);
174  Index n0 = internal::random<Index>(1, outer - (std::max)(j0, j1));
175  if (SparseMatrixType::IsRowMajor)
176  VERIFY_IS_APPROX(m2.innerVectors(j0, n0), refMat2.block(j0, 0, n0, cols));
177  else
178  VERIFY_IS_APPROX(m2.innerVectors(j0, n0), refMat2.block(0, j0, rows, n0));
179  if (SparseMatrixType::IsRowMajor)
180  VERIFY_IS_APPROX(m2.innerVectors(j0, n0) + m2.innerVectors(j1, n0),
181  refMat2.middleRows(j0, n0) + refMat2.middleRows(j1, n0));
182  else
183  VERIFY_IS_APPROX(m2.innerVectors(j0, n0) + m2.innerVectors(j1, n0),
184  refMat2.block(0, j0, rows, n0) + refMat2.block(0, j1, rows, n0));
185 
186  VERIFY_IS_APPROX(m2, refMat2);
187 
188  VERIFY(m2.innerVectors(j0, n0).nonZeros() == m2.transpose().innerVectors(j0, n0).nonZeros());
189 
190  m2.innerVectors(j0, n0) = m2.innerVectors(j0, n0) + m2.innerVectors(j1, n0);
191  if (SparseMatrixType::IsRowMajor)
192  refMat2.middleRows(j0, n0) = (refMat2.middleRows(j0, n0) + refMat2.middleRows(j1, n0)).eval();
193  else
194  refMat2.middleCols(j0, n0) = (refMat2.middleCols(j0, n0) + refMat2.middleCols(j1, n0)).eval();
195 
196  VERIFY_IS_APPROX(m2, refMat2);
197  }
198 
199  // test generic blocks
200  {
201  DenseMatrix refMat2 = DenseMatrix::Zero(rows, cols);
202  SparseMatrixType m2(rows, cols);
203  initSparse<Scalar>(density, refMat2, m2);
204  Index j0 = internal::random<Index>(0, outer - 2);
205  Index j1 = internal::random<Index>(0, outer - 2);
206  Index n0 = internal::random<Index>(1, outer - (std::max)(j0, j1));
207  if (SparseMatrixType::IsRowMajor)
208  VERIFY_IS_APPROX(m2.block(j0, 0, n0, cols), refMat2.block(j0, 0, n0, cols));
209  else
210  VERIFY_IS_APPROX(m2.block(0, j0, rows, n0), refMat2.block(0, j0, rows, n0));
211 
212  if (SparseMatrixType::IsRowMajor)
213  VERIFY_IS_APPROX(m2.block(j0, 0, n0, cols) + m2.block(j1, 0, n0, cols),
214  refMat2.block(j0, 0, n0, cols) + refMat2.block(j1, 0, n0, cols));
215  else
216  VERIFY_IS_APPROX(m2.block(0, j0, rows, n0) + m2.block(0, j1, rows, n0),
217  refMat2.block(0, j0, rows, n0) + refMat2.block(0, j1, rows, n0));
218 
219  Index i = internal::random<Index>(0, m2.outerSize() - 1);
220  if (SparseMatrixType::IsRowMajor) {
221  m2.innerVector(i) = m2.innerVector(i) * s1;
222  refMat2.row(i) = refMat2.row(i) * s1;
223  VERIFY_IS_APPROX(m2, refMat2);
224  } else {
225  m2.innerVector(i) = m2.innerVector(i) * s1;
226  refMat2.col(i) = refMat2.col(i) * s1;
227  VERIFY_IS_APPROX(m2, refMat2);
228  }
229 
230  Index r0 = internal::random<Index>(0, rows - 2);
231  Index c0 = internal::random<Index>(0, cols - 2);
232  Index r1 = internal::random<Index>(1, rows - r0);
233  Index c1 = internal::random<Index>(1, cols - c0);
234 
235  VERIFY_IS_APPROX(DenseVector(m2.col(c0)), refMat2.col(c0));
236  VERIFY_IS_APPROX(m2.col(c0), refMat2.col(c0));
237 
238  VERIFY_IS_APPROX(RowDenseVector(m2.row(r0)), refMat2.row(r0));
239  VERIFY_IS_APPROX(m2.row(r0), refMat2.row(r0));
240 
241  VERIFY_IS_APPROX(m2.block(r0, c0, r1, c1), refMat2.block(r0, c0, r1, c1));
242  VERIFY_IS_APPROX((2 * m2).block(r0, c0, r1, c1), (2 * refMat2).block(r0, c0, r1, c1));
243 
244  if (m2.nonZeros() > 0) {
245  VERIFY_IS_APPROX(m2, refMat2);
246  SparseMatrixType m3(rows, cols);
247  DenseMatrix refMat3(rows, cols);
248  refMat3.setZero();
249  Index n = internal::random<Index>(1, 10);
250  for (Index k = 0; k < n; ++k) {
251  Index o1 = internal::random<Index>(0, outer - 1);
252  Index o2 = internal::random<Index>(0, outer - 1);
253  if (SparseMatrixType::IsRowMajor) {
254  m3.innerVector(o1) = m2.row(o2);
255  refMat3.row(o1) = refMat2.row(o2);
256  } else {
257  m3.innerVector(o1) = m2.col(o2);
258  refMat3.col(o1) = refMat2.col(o2);
259  }
260  if (internal::random<bool>()) m3.makeCompressed();
261  }
262  if (m3.nonZeros() > 0) VERIFY_IS_APPROX(m3, refMat3);
263  }
264  }
265 
266  // Explicit inner iterator.
267  {
268  DenseMatrix refMat2 = DenseMatrix::Zero(rows, cols);
269  SparseMatrixType m2(rows, cols);
270  initSparse<Scalar>(density, refMat2, m2);
271 
272  Index j0 = internal::random<Index>(0, outer - 1);
273  auto v = innervec(m2, j0);
274 
275  typename decltype(v)::InnerIterator block_iterator(v);
276  typename SparseMatrixType::InnerIterator matrix_iterator(m2, j0);
277  while (block_iterator) {
278  VERIFY_IS_EQUAL(block_iterator.index(), matrix_iterator.index());
279  ++block_iterator;
280  ++matrix_iterator;
281  }
282  }
283 }
v
Array< int, Dynamic, 1 > v
Definition: Array_initializer_list_vector_cxx11.cpp:1
DenseMatrix
Matrix< Scalar, Dynamic, Dynamic > DenseMatrix
Definition: BenchSparseUtil.h:23
DenseVector
Matrix< Scalar, Dynamic, 1 > DenseVector
Definition: BenchSparseUtil.h:24
n
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
w
RowVector3d w
Definition: Matrix_resize_int.cpp:3
m2
MatrixType m2(n_dims)
v1
M1<< 1, 2, 3, 4, 5, 6, 7, 8, 9;Map< RowVectorXf > v1(M1.data(), M1.size())
block
m m block(1, 0, 2, 2)<< 4
rows
int rows
Definition: Tutorial_commainit_02.cpp:1
cols
int cols
Definition: Tutorial_commainit_02.cpp:1
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:45
RealScalar
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
Eigen::InnerIterator
An InnerIterator allows to loop over the element of any matrix expression.
Definition: CoreIterators.h:37
Eigen::PlainObjectBase::coeff
EIGEN_DEVICE_FUNC constexpr EIGEN_STRONG_INLINE const Scalar & coeff(Index rowId, Index colId) const
Definition: PlainObjectBase.h:198
Eigen::SparseVector
a sparse vector class
Definition: SparseVector.h:62
real
float real
Definition: datatypes.h:10
min
#define min(a, b)
Definition: datatypes.h:22
max
#define max(a, b)
Definition: datatypes.h:23
VERIFY_IS_APPROX
#define VERIFY_IS_APPROX(a, b)
Definition: integer_types.cpp:13
m
int * m
Definition: level2_cplx_impl.h:294
k
char char char int int * k
Definition: level2_impl.h:374
VERIFY
#define VERIFY(a)
Definition: main.h:362
VERIFY_IS_EQUAL
#define VERIFY_IS_EQUAL(a, b)
Definition: main.h:367
Eigen::numext::is_exactly_zero
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool is_exactly_zero(const X &x)
Definition: Meta.h:592
Eigen::Index
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:83
UniformPSDSelfTest.density
density
Definition: UniformPSDSelfTest.py:19
oomph::PseudoSolidHelper::Zero
double Zero
Definition: pseudosolid_node_update_elements.cc:35
plotPSD.t
t
Definition: plotPSD.py:36
innervec
std::enable_if_t<(T::Flags &RowMajorBit)==RowMajorBit, typename T::RowXpr > innervec(T &A, Index i)
Definition: sparse_block.cpp:14
eval
internal::nested_eval< T, 1 >::type eval(const T &xpr)
Definition: sparse_permutations.cpp:47
j
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References block(), calibrate::c, Eigen::PlainObjectBase< Derived >::coeff(), cols, UniformPSDSelfTest::density, eval(), i, Eigen::InnerIterator< XprType >::index(), innervec(), Eigen::numext::is_exactly_zero(), j, k, m, m2(), max, min, n, UniformPSDSelfTest::r, rows, Eigen::PlainObjectBase< Derived >::setZero(), plotPSD::t, v, v1(), VERIFY, VERIFY_IS_APPROX, VERIFY_IS_EQUAL, w, and oomph::PseudoSolidHelper::Zero.

Referenced by EIGEN_DECLARE_TEST().