From 3d0a3410ad5a23d33c7e445797dba21a4404c674 Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Mon, 6 Jun 2016 17:12:33 +0200
Subject: [PATCH] adapted to changed matrix specs
---
MathLib/LinAlg/MatrixVectorTraits.cpp | 44 ++++++-------------
.../LocalLinearLeastSquaresExtrapolator.h | 17 ++++---
.../GroundwaterFlow/GroundwaterFlowProcess.h | 3 +-
ProcessLib/Process.h | 9 ++--
ProcessLib/TES/TESProcess.cpp | 31 ++++++++-----
Tests/NumLib/ODEs.h | 6 +--
Tests/NumLib/TestExtrapolation.cpp | 8 +++-
7 files changed, 60 insertions(+), 58 deletions(-)
diff --git a/MathLib/LinAlg/MatrixVectorTraits.cpp b/MathLib/LinAlg/MatrixVectorTraits.cpp
index 2ea58fb7851..c2e21c2fe3c 100644
--- a/MathLib/LinAlg/MatrixVectorTraits.cpp
+++ b/MathLib/LinAlg/MatrixVectorTraits.cpp
@@ -34,10 +34,8 @@ std::unique_ptr<Eigen::MatrixXd>
MatrixVectorTraits<Eigen::MatrixXd>::
newInstance(MatrixSpecifications const& spec)
{
- auto const nrows = spec.dof_table ? spec.dof_table->dofSizeWithGhosts() : spec.nrows;
- auto const ncols = spec.dof_table ? nrows : spec.ncols;
-
- return std::unique_ptr<Eigen::MatrixXd>(new Eigen::MatrixXd(nrows, ncols));
+ return std::unique_ptr<Eigen::MatrixXd>(
+ new Eigen::MatrixXd(spec.nrows, spec.ncols));
}
std::unique_ptr<Eigen::VectorXd>
@@ -58,9 +56,7 @@ std::unique_ptr<Eigen::VectorXd>
MatrixVectorTraits<Eigen::VectorXd>::
newInstance(MatrixSpecifications const& spec)
{
- auto const nrows = spec.dof_table ? spec.dof_table->dofSizeWithGhosts() : spec.nrows;
-
- return std::unique_ptr<Eigen::VectorXd>(new Eigen::VectorXd(nrows));
+ return std::unique_ptr<Eigen::VectorXd>(new Eigen::VectorXd(spec.nrows));
}
} // namespace MathLib
@@ -93,22 +89,15 @@ std::unique_ptr<PETScMatrix>
MatrixVectorTraits<PETScMatrix>::
newInstance(MatrixSpecifications const& spec)
{
- auto const nrows =
- spec.dof_table ? spec.dof_table->dofSizeWithoutGhosts() : spec.nrows;
- auto const ncols = spec.dof_table ? nrows : spec.ncols;
+ auto const nrows = spec.nrows;
+ auto const ncols = spec.ncols;
- // TODO I guess it is not hard to make NumLib::computeSparsityPattern()
- // work also for NodePartitionedMesh'es.
- assert(((!spec.sparsity_pattern) || spec.sparsity_pattern->size() == 0u) &&
- "The sparsity pattern is not used with PETSc so I rather crash than"
- " silently accept a precomputed sparsity pattern.");
-
- if (spec.mesh)
+ if (spec.sparsity_pattern)
{
- assert(dynamic_cast<MeshLib::NodePartitionedMesh const*>(spec.mesh));
- auto const& mesh = *static_cast<MeshLib::NodePartitionedMesh const*>(spec.mesh);
+ // Assert that the misuse of the sparsity pattern is consistent.
+ assert(spec.sparsity_pattern->size() == 1);
- auto const max_nonzeroes = mesh.getMaximumNConnectedNodesToNode();
+ auto const max_nonzeroes = spec.sparsity_pattern->front();
PETScMatrixOption mat_opt;
mat_opt.d_nz = max_nonzeroes;
@@ -142,10 +131,9 @@ newInstance(MatrixSpecifications const& spec)
{
auto const is_global_size = false;
- if (spec.dof_table) {
- auto const& dt = *spec.dof_table;
- return std::unique_ptr<PETScVector>(new PETScVector(
- dt.dofSizeWithoutGhosts(), dt.getGhostIndices(), is_global_size));
+ if (spec.ghost_indices != nullptr) {
+ return std::unique_ptr<PETScVector>(
+ new PETScVector(spec.nrows, *spec.ghost_indices, is_global_size));
} else {
return std::unique_ptr<PETScVector>(
new PETScVector(spec.nrows, is_global_size));
@@ -178,9 +166,7 @@ std::unique_ptr<EigenMatrix>
MatrixVectorTraits<EigenMatrix>::
newInstance(MatrixSpecifications const& spec)
{
- auto const nrows = spec.dof_table ? spec.dof_table->dofSizeWithGhosts() : spec.nrows;
-
- auto A = std::unique_ptr<EigenMatrix>(new EigenMatrix(nrows));
+ auto A = std::unique_ptr<EigenMatrix>(new EigenMatrix(spec.nrows));
if (spec.sparsity_pattern)
setMatrixSparsity(*A, *spec.sparsity_pattern);
@@ -206,9 +192,7 @@ std::unique_ptr<EigenVector>
MatrixVectorTraits<EigenVector>::
newInstance(MatrixSpecifications const& spec)
{
- auto const nrows = spec.dof_table ? spec.dof_table->dofSizeWithGhosts() : spec.nrows;
-
- return std::unique_ptr<EigenVector>(new EigenVector(nrows));
+ return std::unique_ptr<EigenVector>(new EigenVector(spec.nrows));
}
} // namespace MathLib
diff --git a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
index b15c46ab429..a2f2e363b96 100644
--- a/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
+++ b/NumLib/Extrapolation/LocalLinearLeastSquaresExtrapolator.h
@@ -49,15 +49,20 @@ public:
* dof_table for one single component variable.
*/
explicit LocalLinearLeastSquaresExtrapolator(
- MathLib::MatrixSpecifications const& matrix_specs)
- : _nodal_values(MathLib::GlobalVectorProvider<GlobalVector>::provider.
- getVector(matrix_specs))
+ MathLib::MatrixSpecifications const& matrix_specs,
+ NumLib::LocalToGlobalIndexMap const& dof_table)
+ : _nodal_values(
+ MathLib::GlobalVectorProvider<GlobalVector>::provider.getVector(
+ matrix_specs))
#ifndef USE_PETSC
- , _residuals(matrix_specs.dof_table->size())
+ ,
+ _residuals(dof_table.size())
#else
- , _residuals(matrix_specs.dof_table->size(), false)
+ ,
+ _residuals(dof_table.size(), false)
#endif
- , _local_to_global(*matrix_specs.dof_table)
+ ,
+ _local_to_global(dof_table)
{}
diff --git a/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h b/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
index 8f6bc793428..6ca68d79db6 100644
--- a/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
+++ b/ProcessLib/GroundwaterFlow/GroundwaterFlowProcess.h
@@ -102,7 +102,8 @@ private:
// TOOD Later on the DOF table can change during the simulation!
_extrapolator.reset(
- new ExtrapolatorImplementation(Base::getMatrixSpecifications()));
+ new ExtrapolatorImplementation(Base::getMatrixSpecifications(),
+ *Base::_local_to_global_index_map));
Base::_secondary_variables.addSecondaryVariable(
"darcy_velocity_x", 1,
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index f543196a2f8..9b97d21e0b8 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -87,10 +87,8 @@ public:
DBUG("Construct dof mappings.");
constructDofTable();
-#ifndef USE_PETSC
DBUG("Compute sparsity pattern");
computeSparsityPattern();
-#endif
initializeConcreteProcess(*_local_to_global_index_map, _mesh,
_integration_order);
@@ -127,10 +125,9 @@ public:
MathLib::MatrixSpecifications getMatrixSpecifications() const override final
{
- return { 0u, 0u,
- &_sparsity_pattern,
- _local_to_global_index_map.get(),
- &_mesh };
+ auto const& l = *_local_to_global_index_map;
+ return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
+ &l.getGhostIndices(), &_sparsity_pattern};
}
void assemble(const double t, GlobalVector const& x,
diff --git a/ProcessLib/TES/TESProcess.cpp b/ProcessLib/TES/TESProcess.cpp
index 48b887a7060..b1bf3d20c4d 100644
--- a/ProcessLib/TES/TESProcess.cpp
+++ b/ProcessLib/TES/TESProcess.cpp
@@ -192,10 +192,15 @@ void TESProcess<GlobalSetup>::initializeConcreteProcess(
// by location order is needed for output
NumLib::ComponentOrder::BY_LOCATION));
- _extrapolator.reset(
- new ExtrapolatorImplementation(MathLib::MatrixSpecifications(
- 0u, 0u, nullptr, _local_to_global_index_map_single_component.get(),
- &mesh)));
+ {
+ auto const& l = *_local_to_global_index_map_single_component;
+ _extrapolator.reset(new ExtrapolatorImplementation(
+ MathLib::MatrixSpecifications(l.dofSizeWithoutGhosts(),
+ l.dofSizeWithoutGhosts(),
+ &l.getGhostIndices(),
+ nullptr),
+ l));
+ }
// secondary variables
auto add2nd = [&](std::string const& var_name, unsigned const n_components,
@@ -343,9 +348,11 @@ TESProcess<GlobalSetup>::computeVapourPartialPressure(
{
assert(&dof_table == this->_local_to_global_index_map.get());
+ auto const& dof_table_single = *_local_to_global_index_map_single_component;
result_cache = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(
- {0, 0, nullptr, _local_to_global_index_map_single_component.get(),
- nullptr});
+ {dof_table_single.dofSizeWithoutGhosts(),
+ dof_table_single.dofSizeWithoutGhosts(),
+ &dof_table_single.getGhostIndices(), nullptr});
GlobalIndexType nnodes = this->_mesh.getNNodes();
@@ -375,9 +382,11 @@ TESProcess<GlobalSetup>::computeRelativeHumidity(
{
assert(&dof_table == this->_local_to_global_index_map.get());
+ auto const& dof_table_single = *_local_to_global_index_map_single_component;
result_cache = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(
- {0, 0, nullptr, _local_to_global_index_map_single_component.get(),
- nullptr});
+ {dof_table_single.dofSizeWithoutGhosts(),
+ dof_table_single.dofSizeWithoutGhosts(),
+ &dof_table_single.getGhostIndices(), nullptr});
GlobalIndexType nnodes = this->_mesh.getNNodes();
@@ -412,9 +421,11 @@ TESProcess<GlobalSetup>::computeEquilibriumLoading(
{
assert(&dof_table == this->_local_to_global_index_map.get());
+ auto const& dof_table_single = *_local_to_global_index_map_single_component;
result_cache = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(
- {0, 0, nullptr, _local_to_global_index_map_single_component.get(),
- nullptr});
+ {dof_table_single.dofSizeWithoutGhosts(),
+ dof_table_single.dofSizeWithoutGhosts(),
+ &dof_table_single.getGhostIndices(), nullptr});
GlobalIndexType nnodes = this->_mesh.getNNodes();
diff --git a/Tests/NumLib/ODEs.h b/Tests/NumLib/ODEs.h
index b693bf6efdb..e864206aab7 100644
--- a/Tests/NumLib/ODEs.h
+++ b/Tests/NumLib/ODEs.h
@@ -55,7 +55,7 @@ public:
MathLib::MatrixSpecifications getMatrixSpecifications() const override
{
- return { N, N, nullptr, nullptr, nullptr };
+ return { N, N, nullptr, nullptr };
}
bool isLinear() const override
@@ -133,7 +133,7 @@ public:
MathLib::MatrixSpecifications getMatrixSpecifications() const override
{
- return { N, N, nullptr, nullptr, nullptr };
+ return { N, N, nullptr, nullptr };
}
bool isLinear() const override
@@ -271,7 +271,7 @@ public:
MathLib::MatrixSpecifications getMatrixSpecifications() const override
{
- return { N, N, nullptr, nullptr, nullptr };
+ return { N, N, nullptr, nullptr };
}
bool isLinear() const override
diff --git a/Tests/NumLib/TestExtrapolation.cpp b/Tests/NumLib/TestExtrapolation.cpp
index c7afac39d8f..62262b063c2 100644
--- a/Tests/NumLib/TestExtrapolation.cpp
+++ b/Tests/NumLib/TestExtrapolation.cpp
@@ -182,7 +182,11 @@ public:
// Passing _dof_table works, because this process has only one variable
// and the variable has exactly one component.
_extrapolator.reset(new ExtrapolatorImplementation(
- MathLib::MatrixSpecifications{0u, 0u, nullptr, _dof_table.get(), nullptr}));
+ MathLib::MatrixSpecifications{_dof_table->dofSizeWithoutGhosts(),
+ _dof_table->dofSizeWithoutGhosts(),
+ &_dof_table->getGhostIndices(),
+ nullptr},
+ *_dof_table));
createAssemblers(mesh);
}
@@ -333,7 +337,7 @@ TEST(NumLib, DISABLED_Extrapolation)
TestProcess<GlobalSetup> pcs(*mesh, integration_order);
// generate random nodal values
- MathLib::MatrixSpecifications spec{nnodes, nnodes, nullptr, nullptr, nullptr};
+ MathLib::MatrixSpecifications spec{nnodes, nnodes, nullptr, nullptr};
auto x = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(spec);
fillVectorRandomly(*x);
--
GitLab