From c383fc2f40f57dc111329ddeaa726f29446a7714 Mon Sep 17 00:00:00 2001
From: Patrick Jaap <jaap@wias-berlin.de>
Date: Tue, 2 May 2023 10:36:25 +0200
Subject: [PATCH] Prepare UMFPACK for MultiTypeBlockMatrix
---
CHANGELOG.md | 4 ++++
dune/istl/umfpack.hh | 21 +++++++++++++++++++++
2 files changed, 25 insertions(+)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 1fc0d4a49..d037dec9a 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -5,6 +5,10 @@ SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
# Master (will become release 2.10)
+- `UMFPACK` is extended to arbitrary blocked matrix structures. This includes `MultiTypeBlockMatrix`. The external interface is unchanged.
+ However, internally the `BCCSMatrixInitializer` is replaced by direct calls of `flatMatrixForEach` similar to `Cholmod`. This requires a
+ compatible vector field of "ignored" degrees of freedom. The method `setSubMatrix` with a top-level index set is preserved for compatibility.
+
- The internal storage in `MatrixIndexSet` was changed from `std::set` to a sorted `std::vector`
(with a `std::set` fallback for very dense rows) to improve performance. The stored index
type was changed from `std::size_t` to `uint32_t` to reduce memory consumption and improve
diff --git a/dune/istl/umfpack.hh b/dune/istl/umfpack.hh
index 6f2e8771b..63fa05dcd 100644
--- a/dune/istl/umfpack.hh
+++ b/dune/istl/umfpack.hh
@@ -18,6 +18,8 @@
#include<dune/istl/bccsmatrixinitializer.hh>
#include<dune/istl/bcrsmatrix.hh>
#include<dune/istl/foreach.hh>
+#include<dune/istl/multitypeblockmatrix.hh>
+#include<dune/istl/multitypeblockvector.hh>
#include<dune/istl/solvers.hh>
#include<dune/istl/solvertype.hh>
#include <dune/istl/solverfactory.hh>
@@ -195,6 +197,25 @@ namespace Dune {
using range_type = BlockVector<T, A>;
};
+ // to make the `UMFPackVectorChooser` work with `MultiTypeBlockMatrix`, we need to add an intermediate step for the rows, which are typically `MultiTypeBlockVector`
+ template<typename FirstBlock, typename... Blocks>
+ struct UMFPackVectorChooser<MultiTypeBlockVector<FirstBlock, Blocks...> >
+ {
+ /** @brief The type of the domain of the solver */
+ using domain_type = MultiTypeBlockVector<typename UMFPackVectorChooser<FirstBlock>::domain_type,typename UMFPackVectorChooser<Blocks>::domain_type... >;
+ /** @brief The type of the range of the solver */
+ using range_type = typename UMFPackVectorChooser<FirstBlock>::range_type;
+ };
+
+ // specialization for `MultiTypeBlockMatrix` with `MultiTypeBlockVector` rows
+ template<typename FirstRow, typename... Rows>
+ struct UMFPackVectorChooser<MultiTypeBlockMatrix<FirstRow, Rows...> >
+ {
+ /** @brief The type of the domain of the solver */
+ using domain_type = typename UMFPackVectorChooser<FirstRow>::domain_type;
+ /** @brief The type of the range of the solver */
+ using range_type = MultiTypeBlockVector< typename UMFPackVectorChooser<FirstRow>::range_type, typename UMFPackVectorChooser<Rows>::range_type... >;
+ };
// dummy class to represent no BitVector
struct NoBitVector
--
GitLab