From f1e532ba97cb500a0084ed962ce777ee43d81db4 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Wed, 27 Jul 2016 19:30:20 +0200
Subject: [PATCH] Tests: Use a sparse matrix
Otherwise, rounding errors dominate once the problem becomes larger
---
.../matrix-vector/test/triangularsolvetest.cc | 35 +++++++++++++------
1 file changed, 24 insertions(+), 11 deletions(-)
diff --git a/dune/matrix-vector/test/triangularsolvetest.cc b/dune/matrix-vector/test/triangularsolvetest.cc
index 28dd685..bca54c4 100644
--- a/dune/matrix-vector/test/triangularsolvetest.cc
+++ b/dune/matrix-vector/test/triangularsolvetest.cc
@@ -16,34 +16,47 @@
#include "common.hh"
double const tol = 1e-10;
+double const reg = 1e-6;
template <int n, bool lower>
bool test() {
bool passed = true;
+ std::random_device randomDevice;
+ std::mt19937 generator(randomDevice());
+
Dune::BCRSMatrix<Dune::FieldMatrix<double, 1, 1>> M;
Dune::MatrixIndexSet indices(n, n);
- for (size_t i = 0; i < n; ++i)
- for (size_t j = lower ? 0 : i; j < (lower ? i + 1 : n); ++j)
- indices.add(i, j);
+ for (size_t i = 0; i < n; ++i) {
+ size_t initialJ = lower ? 0 : i;
+ size_t finalJ = lower ? i + 1 : n;
+ for (size_t j = initialJ; j < finalJ; ++j) {
+ std::uniform_int_distribution<int> integerDistribution(1, finalJ - initialJ);
+ auto aroundTen = [&]() {
+ auto x = integerDistribution(generator);
+ return x <= 10;
+ };
+ // Add a diagonal entry and approximately ten non-zero entries
+ if (i == j or aroundTen())
+ indices.add(i, j);
+ }
+ }
indices.exportIdx(M);
- std::random_device randomDevice;
- std::mt19937 generator(randomDevice());
- std::uniform_real_distribution<> distribution(0.1, 0.4);
-
+ std::uniform_real_distribution<> lowDistribution(reg, 0.4);
for (auto it = M.begin(); it != M.end(); ++it) {
size_t const i = it.index();
for (auto cIt = it->begin(); cIt != it->end(); ++cIt) {
size_t const j = cIt.index();
- *cIt = distribution(generator);
+ *cIt = lowDistribution(generator);
if (i == j)
*cIt += 0.6;
}
}
+ std::uniform_real_distribution<> unconstrainedDistribution(reg, 0.4);
Dune::BlockVector<Dune::FieldVector<double, 1>> sol(n);
for (auto &x : sol)
- x = distribution(generator);
+ x = unconstrainedDistribution(generator);
Dune::BlockVector<Dune::FieldVector<double, 1>> b(n);
M.mv(sol, b);
@@ -72,7 +85,7 @@ bool test() {
int main() {
bool passed = true;
- passed &= test<1000, false>();
- passed &= test<1000, true>();
+ passed &= test<5000, false>();
+ passed &= test<5000, true>();
return passed ? 0 : 1;
}
--
GitLab