From 84d345314c9396c2905f1577fd5982a658557e8a Mon Sep 17 00:00:00 2001
From: podlesny <podlesny@zedat.fu-berlin.de>
Date: Thu, 8 Mar 2018 15:59:51 +0100
Subject: [PATCH] .
---
src/explicitgrid.hh | 6 ++++
src/multi-body-problem.cc | 19 +++++++++--
src/spatial-solving/fixedpointiterator.cc | 33 +++++++++----------
src/spatial-solving/fixedpointiterator.hh | 13 ++++++--
.../fixedpointiterator_tmpl.cc | 4 ++-
src/spatial-solving/solverfactory.cc | 16 ++++-----
src/spatial-solving/solverfactory.hh | 6 ++--
7 files changed, 64 insertions(+), 33 deletions(-)
diff --git a/src/explicitgrid.hh b/src/explicitgrid.hh
index b538e24d..685d2a49 100644
--- a/src/explicitgrid.hh
+++ b/src/explicitgrid.hh
@@ -2,7 +2,13 @@
#define SRC_EXPLICITGRID_HH
#include "gridselector.hh"
+#include "explicitvectors.hh"
+
+#include <dune/contact/common/deformedcontinuacomplex.hh>
using LeafGridView = Grid::LeafGridView;
using LevelGridView = Grid::LevelGridView;
+
+using DeformedGridComplex = typename Dune::Contact::DeformedContinuaComplex<Grid, Vector>;
+using DeformedGrid = DeformedGridComplex::DeformedGridType;
#endif
diff --git a/src/multi-body-problem.cc b/src/multi-body-problem.cc
index 83d79227..9f2b6c67 100644
--- a/src/multi-body-problem.cc
+++ b/src/multi-body-problem.cc
@@ -65,11 +65,11 @@
#include "multi-body-problem-data/myglobalfrictiondata.hh"
#include "multi-body-problem-data/mygrids.hh"
#include "multi-body-problem-data/weakpatch.hh"
-#include "spatial-solving/solverfactory.hh"
-#include "time-stepping/adaptivetimestepper.hh"
+//#include "spatial-solving/solverfactory.hh"
+//#include "time-stepping/adaptivetimestepper.hh"
#include "time-stepping/rate.hh"
#include "time-stepping/state.hh"
-#include "time-stepping/updaters.hh"
+//#include "time-stepping/updaters.hh"
#include "vtk.hh"
// for getcwd
@@ -194,8 +194,21 @@ int main(int argc, char *argv[]) {
for (size_t i=0; i<grids.size(); i++) {
deformedGridComplex.addGrid(*grids[i], "body" + std::to_string(i));
}
+
const auto deformedGrids = deformedGridComplex.gridVector();
+ /*
+ // test deformed grids
+ for (size_t i=0; i<deformedGrids.size(); i++) {
+ Vector def(deformedGrids[i]->size(dims));
+ def = 1;
+ deformedGridComplex.setDeformation(def, i);
+
+ writeToVTK(deformedGrids[i]->leafGridView(), "", "deformedGrid" + std::to_string(i));
+ }
+ // test deformed grids
+ */
+
std::vector<std::shared_ptr<DefLeafGridView>> leafViews(deformedGrids.size());
std::vector<std::shared_ptr<DefLevelGridView>> levelViews(deformedGrids.size());
std::vector<int> leafVertexCounts(deformedGrids.size());
diff --git a/src/spatial-solving/fixedpointiterator.cc b/src/spatial-solving/fixedpointiterator.cc
index 02322caf..54a251f0 100644
--- a/src/spatial-solving/fixedpointiterator.cc
+++ b/src/spatial-solving/fixedpointiterator.cc
@@ -32,8 +32,8 @@ void FixedPointIterationCounter::operator+=(
multigridIterations += other.multigridIterations;
}
-template <class DeformedGrid, class Factory, class Updaters, class ErrorNorm>
-FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::FixedPointIterator(
+template <class Factory, class Updaters, class ErrorNorm>
+FixedPointIterator<Factory, Updaters, ErrorNorm>::FixedPointIterator(
const Dune::Contact::NBodyAssembler<DeformedGrid, Vector>& nBodyAssembler,
Factory &factory, Dune::ParameterTree const &parset,
std::shared_ptr<Nonlinearity> globalFriction, ErrorNorm const &errorNorm)
@@ -49,9 +49,9 @@ FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::FixedPointIterat
verbosity_(parset.get<Solver::VerbosityMode>("v.solver.verbosity")),
errorNorm_(errorNorm) {}
-template <class DeformedGrid,class Factory, class Updaters, class ErrorNorm>
+template <class Factory, class Updaters, class ErrorNorm>
FixedPointIterationCounter
-FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::run(
+FixedPointIterator<Factory, Updaters, ErrorNorm>::run(
Updaters updaters, Matrix const &velocityMatrix, Vector const &velocityRHS,
Vector &velocityIterate) {
EnergyNorm<Matrix, Vector> energyNorm(velocityMatrix);
@@ -66,7 +66,7 @@ FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::run(
for (fixedPointIteration = 0; fixedPointIteration < fixedPointMaxIterations_;
++fixedPointIteration) {
// solve a velocity problem
- globalFriction_->updateAlpha(alpha);
+ // globalFriction_->updateAlpha(alpha);
ConvexProblem convexProblem(1.0, velocityMatrix, *globalFriction_,
velocityRHS, velocityIterate);
BlockProblem velocityProblem(parset_, convexProblem);
@@ -82,7 +82,7 @@ FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::run(
for (size_t i=0; i<v_m.size(); i++) {
v_m[i] *= 1.0 - lambda_;
- Arithmetic::addProduct(v_m[i], lambda_, velocityIterate[i]);
+ //Arithmetic::addProduct(v_m[i], lambda_, velocityIterate[i]);
}
// compute relative velocities on contact boundaries
@@ -91,14 +91,14 @@ FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::run(
// solve a state problem
updaters.state_->solve(v_m);
ScalarVector newAlpha;
- updaters.state_->extractAlpha(newAlpha);
+ /* updaters.state_->extractAlpha(newAlpha);
if (lambda_ < 1e-12 or
errorNorm_.diff(alpha, newAlpha) < fixedPointTolerance_) {
fixedPointIteration++;
break;
}
- alpha = newAlpha;
+ alpha = newAlpha;*/
}
if (fixedPointIteration == fixedPointMaxIterations_)
DUNE_THROW(Dune::Exception, "FPI failed to converge");
@@ -119,8 +119,10 @@ std::ostream &operator<<(std::ostream &stream,
<< ")";
}
-template <class DeformedGrid, class Factory, class Updaters, class ErrorNorm>
-void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVelocities(std::vector<Vector>& v_m) const {
+template <class Factory, class Updaters, class ErrorNorm>
+void FixedPointIterator<Factory, Updaters, ErrorNorm>::relativeVelocities(std::vector<Vector>& v_m) const {
+ using field_type = typename Factory::Matrix::field_type;
+
// adaptation of DualMortarCoupling::setup()
const size_t dim = DeformedGrid::dimension;
@@ -134,7 +136,7 @@ void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVel
using DualCache = Dune::Contact::DualBasisAdapter<GridView, field_type>;
std::unique_ptr<DualCache> dualCache;
dualCache = std::make_unique< Dune::Contact::DualBasisAdapterGlobal<GridView, field_type> >();
-
+/*
// define FE grid functions
std::vector<BasisGridFunction<> > gridFunctions(nBodyAssembler_.nGrids());
for (size_t i=0; i<gridFunctions.size(); i++) {
@@ -245,9 +247,6 @@ void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVel
}
-
-
-
///////////////////////////////////
// reducing nonmortar boundary
/////////////////////////////////
@@ -278,7 +277,7 @@ void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVel
//writeBoundary(boundaryWithMapping,debugPath_ + "relevantNonmortar");
- /** \todo replace by all fine grid segments which are totally covered by the RemoteIntersections. */
+ // \todo replace by all fine grid segments which are totally covered by the RemoteIntersections.
//for (const auto& rIs : intersections(glue))
// boundaryWithMapping.addFace(rIs.inside(),rIs.indexInInside());
@@ -303,7 +302,7 @@ void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVel
// Get the occupation structure for the mortar matrix
// ///////////////////////////////////////////////////////////
- /** \todo Also restrict mortar indices and don't use the whole grid level. */
+ // todo Also restrict mortar indices and don't use the whole grid level.
Dune::MatrixIndexSet mortarIndices(nonmortarBoundary_.numVertices(), grid1_->size(dim));
// Create mapping from the global set of block dofs to the ones on the contact boundary
@@ -382,7 +381,7 @@ void FixedPointIterator<DeformedGrid, Factory, Updaters, ErrorNorm>::relativeVel
for(size_t rowIdx=0; rowIdx<weakObstacle_.size(); rowIdx++)
weakObstacle_[rowIdx] *= D[rowIdx][rowIdx][0][0];
- gridGlueBackend_->clear();
+ gridGlueBackend_->clear(); */
}
#include "fixedpointiterator_tmpl.cc"
diff --git a/src/spatial-solving/fixedpointiterator.hh b/src/spatial-solving/fixedpointiterator.hh
index ba31c81e..873b9b42 100644
--- a/src/spatial-solving/fixedpointiterator.hh
+++ b/src/spatial-solving/fixedpointiterator.hh
@@ -8,6 +8,8 @@
#include <dune/solvers/norms/norm.hh>
#include <dune/solvers/solvers/solver.hh>
+#include <dune/contact/assemblers/nbodyassembler.hh>
+
struct FixedPointIterationCounter {
size_t iterations = 0;
size_t multigridIterations = 0;
@@ -27,10 +29,16 @@ class FixedPointIterator {
using BlockProblem = typename Factory::BlockProblem;
using Nonlinearity = typename ConvexProblem::NonlinearityType;
+ using DeformedGrid = typename Factory::DeformedGrid;
+
+private:
+ void relativeVelocities(std::vector<Vector>& v_m) const;
+
public:
- FixedPointIterator(Factory &factory, Dune::ParameterTree const &parset,
+ FixedPointIterator(const Dune::Contact::NBodyAssembler<DeformedGrid, Vector>& nBodyAssembler,
+ Factory &factory, const Dune::ParameterTree& parset,
std::shared_ptr<Nonlinearity> globalFriction,
- ErrorNorm const &errorNorm_);
+ const ErrorNorm& errorNorm);
FixedPointIterationCounter run(Updaters updaters,
Matrix const &velocityMatrix,
@@ -38,6 +46,7 @@ class FixedPointIterator {
Vector &velocityIterate);
private:
+ const Dune::Contact::NBodyAssembler<DeformedGrid, Vector>& nBodyAssembler_;
std::shared_ptr<typename Factory::Step> step_;
Dune::ParameterTree const &parset_;
std::shared_ptr<Nonlinearity> globalFriction_;
diff --git a/src/spatial-solving/fixedpointiterator_tmpl.cc b/src/spatial-solving/fixedpointiterator_tmpl.cc
index a50b93d4..5f4ad7c1 100644
--- a/src/spatial-solving/fixedpointiterator_tmpl.cc
+++ b/src/spatial-solving/fixedpointiterator_tmpl.cc
@@ -10,6 +10,8 @@
#include <dune/solvers/norms/energynorm.hh>
#include <dune/tnnmg/problem-classes/convexproblem.hh>
+#include <dune/contact/common/deformedcontinuacomplex.hh>
+
#include <dune/tectonic/globalfriction.hh>
#include <dune/tectonic/myblockproblem.hh>
@@ -22,7 +24,7 @@ using Function = Dune::VirtualFunction<double, double>;
using Factory = SolverFactory<
MY_DIM,
MyBlockProblem<ConvexProblem<GlobalFriction<Matrix, Vector>, Matrix>>,
- Grid>;
+ DeformedGrid>;
using MyStateUpdater = StateUpdater<ScalarVector, Vector>;
using MyRateUpdater = RateUpdater<Vector, Matrix, Function, MY_DIM>;
using MyUpdaters = Updaters<MyRateUpdater, MyStateUpdater>;
diff --git a/src/spatial-solving/solverfactory.cc b/src/spatial-solving/solverfactory.cc
index 3242724a..88d4d0c0 100644
--- a/src/spatial-solving/solverfactory.cc
+++ b/src/spatial-solving/solverfactory.cc
@@ -11,9 +11,9 @@
#include "solverfactory.hh"
-template <size_t dim, class BlockProblem, class Grid>
-SolverFactory<dim, BlockProblem, Grid>::SolverFactory(
- Dune::ParameterTree const &parset, Grid const &grid,
+template <size_t dim, class BlockProblem, class DeformedGrid>
+SolverFactory<dim, BlockProblem, DeformedGrid>::SolverFactory(
+ Dune::ParameterTree const &parset, const DeformedGrid& grid,
Dune::BitSetVector<dim> const &ignoreNodes)
: baseEnergyNorm(linearBaseSolverStep),
linearBaseSolver(&linearBaseSolverStep,
@@ -33,7 +33,7 @@ SolverFactory<dim, BlockProblem, Grid>::SolverFactory(
// transfer operators
for (auto &&x : transferOperators)
x = new CompressedMultigridTransfer<Vector>;
- TransferOperatorAssembler<Grid>(grid)
+ TransferOperatorAssembler<DeformedGrid>(grid)
.assembleOperatorPointerHierarchy(transferOperators);
linearIterationStep.setTransferOperators(transferOperators);
@@ -44,14 +44,14 @@ SolverFactory<dim, BlockProblem, Grid>::SolverFactory(
multigridStep->ignoreNodes_ = &ignoreNodes;
}
-template <size_t dim, class BlockProblem, class Grid>
-SolverFactory<dim, BlockProblem, Grid>::~SolverFactory() {
+template <size_t dim, class BlockProblem, class DeformedGrid>
+SolverFactory<dim, BlockProblem, DeformedGrid>::~SolverFactory() {
for (auto &&x : transferOperators)
delete x;
}
-template <size_t dim, class BlockProblem, class Grid>
-auto SolverFactory<dim, BlockProblem, Grid>::getStep()
+template <size_t dim, class BlockProblem, class DeformedGrid>
+auto SolverFactory<dim, BlockProblem, DeformedGrid>::getStep()
-> std::shared_ptr<Step> {
return multigridStep;
}
diff --git a/src/spatial-solving/solverfactory.hh b/src/spatial-solving/solverfactory.hh
index 8f05a572..813d2d03 100644
--- a/src/spatial-solving/solverfactory.hh
+++ b/src/spatial-solving/solverfactory.hh
@@ -12,7 +12,7 @@
#include <dune/tnnmg/iterationsteps/genericnonlineargs.hh>
#include <dune/tnnmg/iterationsteps/tnnmgstep.hh>
-template <size_t dim, class BlockProblemTEMPLATE, class Grid>
+template <size_t dim, class BlockProblemTEMPLATE, class DeformedGridTEMPLATE>
class SolverFactory {
public:
using BlockProblem = BlockProblemTEMPLATE;
@@ -20,6 +20,8 @@ class SolverFactory {
using Vector = typename BlockProblem::VectorType;
using Matrix = typename BlockProblem::MatrixType;
+ using DeformedGrid = DeformedGridTEMPLATE;
+
private:
using NonlinearSmoother = GenericNonlinearGS<BlockProblem>;
@@ -27,7 +29,7 @@ class SolverFactory {
using Step =
TruncatedNonsmoothNewtonMultigrid<BlockProblem, NonlinearSmoother>;
- SolverFactory(Dune::ParameterTree const &parset, Grid const &grid,
+ SolverFactory(Dune::ParameterTree const &parset, const DeformedGrid& grid,
Dune::BitSetVector<dim> const &ignoreNodes);
~SolverFactory();
--
GitLab