From 4b9841423b9403b34b1ecd22ac657a7f70d69976 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Thu, 28 Jul 2016 01:06:58 +0200
Subject: [PATCH] Tests: More general test; more helpful output
---
.../matrix-vector/test/triangularsolvetest.cc | 27 ++++++++++++++-----
1 file changed, 21 insertions(+), 6 deletions(-)
diff --git a/dune/matrix-vector/test/triangularsolvetest.cc b/dune/matrix-vector/test/triangularsolvetest.cc
index bca54c4..fec0232 100644
--- a/dune/matrix-vector/test/triangularsolvetest.cc
+++ b/dune/matrix-vector/test/triangularsolvetest.cc
@@ -11,6 +11,7 @@
#include <dune/istl/bvector.hh>
#include <dune/istl/matrixindexset.hh>
+#include "../transpose.hh"
#include "../triangularsolve.hh"
#include "common.hh"
@@ -42,18 +43,18 @@ bool test() {
}
indices.exportIdx(M);
- std::uniform_real_distribution<> lowDistribution(reg, 0.4);
+ std::uniform_real_distribution<> unconstrainedDistribution(-1, 1);
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 = lowDistribution(generator);
+ *cIt = unconstrainedDistribution(generator);
+ // Make sure the diagonal is bounded away from zero
if (i == j)
- *cIt += 0.6;
+ *cIt += (*cIt > 0 ? 1 : -1) + reg;
}
}
- std::uniform_real_distribution<> unconstrainedDistribution(reg, 0.4);
Dune::BlockVector<Dune::FieldVector<double, 1>> sol(n);
for (auto &x : sol)
x = unconstrainedDistribution(generator);
@@ -67,7 +68,7 @@ bool test() {
Dune::BlockVector<Dune::FieldVector<double, 1>> y(n);
M.mv(x, y);
auto diff = diffDune(b, y);
- std::cout << "Difference: " << diff << std::endl;
+ std::cout << "|x - A^{-1}b| = " << diff << std::endl;
passed &= diff < tol;
}
{
@@ -77,7 +78,19 @@ bool test() {
Dune::BlockVector<Dune::FieldVector<double, 1>> y(n);
M.mtv(x, y);
auto diff = diffDune(b, y);
- std::cout << "Difference: " << diff << std::endl;
+ std::cout << "|x - A^{-T}b| = " << diff << " (manual transpose)" << std::endl;
+ passed &= diff < tol;
+ }
+ {
+ Dune::BCRSMatrix<Dune::FieldMatrix<double, 1, 1>> Mt;
+ Dune::MatrixVector::transpose(M, Mt);
+ auto x = lower
+ ? Dune::MatrixVector::upperTriangularSolve(Mt, b, ignore)
+ : Dune::MatrixVector::lowerTriangularSolve(Mt, b, ignore);
+ Dune::BlockVector<Dune::FieldVector<double, 1>> y(n);
+ M.mtv(x, y);
+ auto diff = diffDune(b, y);
+ std::cout << "|x - A^{-T}b| = " << diff << std::endl;
passed &= diff < tol;
}
return passed;
@@ -85,7 +98,9 @@ bool test() {
int main() {
bool passed = true;
+ std::cout << "Testing random upper triangular matrix" << std::endl;
passed &= test<5000, false>();
+ std::cout << "Testing lower triangular matrix" << std::endl;
passed &= test<5000, true>();
return passed ? 0 : 1;
}
--
GitLab