From 1ef5b30f47b7086d538f45e6f8438abfbf96e761 Mon Sep 17 00:00:00 2001
From: podlesny <podlesny@mi.fu-berlin.de>
Date: Fri, 26 Jul 2019 17:20:34 +0200
Subject: [PATCH] .
---
src/multi-body-problem-2D.cfg | 2 +-
src/multi-body-problem.cfg | 3 ++-
src/solverfactorytest.cc | 40 +++++++++++++++++-----------
src/spatial-solving/solverfactory.cc | 4 ---
src/spatial-solving/solverfactory.hh | 2 +-
5 files changed, 29 insertions(+), 22 deletions(-)
diff --git a/src/multi-body-problem-2D.cfg b/src/multi-body-problem-2D.cfg
index a308c80b..647e99cb 100644
--- a/src/multi-body-problem-2D.cfg
+++ b/src/multi-body-problem-2D.cfg
@@ -1,6 +1,6 @@
# -*- mode:conf -*-
[boundary.friction]
-smallestDiameter = 0.01 # 2e-3 [m]
+smallestDiameter = 0.05 # 2e-3 [m]
[timeSteps]
refinementTolerance = 1e-5 # 1e-5
diff --git a/src/multi-body-problem.cfg b/src/multi-body-problem.cfg
index e8a3e97b..8df2a20d 100644
--- a/src/multi-body-problem.cfg
+++ b/src/multi-body-problem.cfg
@@ -47,6 +47,7 @@ relativeTau = 1e-4 # 1e-6
[timeSteps]
scheme = newmark
+timeSteps = 5
[u0.solver]
maximumIterations = 100
@@ -76,7 +77,7 @@ maximumIterations = 10000
verbosity = quiet
[solver.tnnmg.main]
-pre = 1
+pre = 10
multi = 5 # number of multigrid steps
post = 0
diff --git a/src/solverfactorytest.cc b/src/solverfactorytest.cc
index 1d0f983b..9ae2798e 100644
--- a/src/solverfactorytest.cc
+++ b/src/solverfactorytest.cc
@@ -57,8 +57,6 @@
const int dim = MY_DIM;
-const int timeSteps = 5;
-
Dune::ParameterTree parset;
size_t bodyCount;
@@ -71,6 +69,7 @@ std::vector<ScalarVector> weightedNormalStress;
double relativeTime;
double relativeTau;
size_t timeStep = 0;
+size_t timeSteps = 0;
const std::string path = "";
const std::string outputFile = "solverfactorytest.log";
@@ -78,7 +77,7 @@ const std::string outputFile = "solverfactorytest.log";
template <class ContactNetwork>
void solveProblem(const ContactNetwork& contactNetwork,
const Matrix& mat, const Vector& rhs, Vector& x,
- const BitVector& ignore, const Vector& lower, const Vector& upper) {
+ const BitVector& ignore, const Vector& lower, const Vector& upper, bool initial = false) {
using Solver = typename Dune::Solvers::LoopSolver<Vector, BitVector>;
using Norm = EnergyNorm<Matrix, Vector>;
@@ -96,6 +95,8 @@ void solveProblem(const ContactNetwork& contactNetwork,
using ContactFunctional = Dune::TNNMG::BoxConstrainedQuadraticFunctional<Matrix&, Vector&, Vector&, Vector&, double>;
auto I = ContactFunctional(mat, rhs, lower, upper);
+ std::cout << "Energy start iterate: " << I(x) << std::endl;
+
using LocalSolver = Dune::TNNMG::ScalarObstacleSolver;
auto localSolver = Dune::TNNMG::gaussSeidelLocalSolver(LocalSolver());
@@ -176,10 +177,16 @@ void solveProblem(const ContactNetwork& contactNetwork,
refSolver.solve();
std::cout << correctionNorms.size() << std::endl;
-
+ if (initial) {
+ x = refX;
+ return;
+ }
// set up solver factory solver
// set up functional
+ /*auto& globalFriction = contactNetwork.globalFriction();
+ using MyFunctional = Functional<Matrix&, Vector&, std::decay_t<decltype(globalFriction)>&, Vector&, Vector&, typename Matrix::field_type>;
+ MyFunctional J(mat, rhs, globalFriction, lower, upper);*/
using MyFunctional = Functional<Matrix&, Vector&, ZeroNonlinearity&, Vector&, Vector&, typename Matrix::field_type>;
MyFunctional J(mat, rhs, ZeroNonlinearity(), lower, upper);
@@ -216,7 +223,7 @@ void solveProblem(const ContactNetwork& contactNetwork,
diff -= refX;
std::cout << "Solver error in energy norm: " << norm(diff) << std::endl;
- x = refX;
+ std::cout << "Energy end iterate: " << J(x) << std::endl;
}
@@ -274,7 +281,7 @@ void setupInitialConditions(const ContactNetwork& contactNetwork) {
print(lower, "lower");
print(upper, "upper");*/
- solveProblem(contactNetwork, bilinearForm, totalRhs, totalX, _dirichletNodes, lower, upper);
+ solveProblem(contactNetwork, bilinearForm, totalRhs, totalX, _dirichletNodes, lower, upper, true);
nBodyAssembler.postprocess(totalX, _x);
};
@@ -434,17 +441,17 @@ void run(Updaters& updaters, ContactNetwork& contactNetwork,
}
// compute velocity obstacles
- Vector vLower, vUpper;
+ /*Vector vLower, vUpper;
std::vector<Vector> u0, v0;
updaters.rate_->extractOldVelocity(v0);
updaters.rate_->extractOldDisplacement(u0);
Vector totalu0, totalv0;
- nBodyAssembler.concatenateVectors(u0, totalu0);
- nBodyAssembler.concatenateVectors(v0, totalv0);
+ nBodyAssembler.nodalToTransformed(u0, totalu0);
+ nBodyAssembler.nodalToTransformed(v0, totalv0);
- //updaters.rate_->velocityObstacles(totalu0, lower, totalv0, vLower);
- //updaters.rate_->velocityObstacles(totalu0, upper, totalv0, vUpper);
+ updaters.rate_->velocityObstacles(totalu0, lower, totalv0, vLower);
+ updaters.rate_->velocityObstacles(totalu0, upper, totalv0, vUpper);*/
const auto& errorNorms = contactNetwork.stateEnergyNorms();
@@ -459,22 +466,24 @@ void run(Updaters& updaters, ContactNetwork& contactNetwork,
size_t fixedPointIteration;
size_t multigridIterations = 0;
+
std::vector<ScalarVector> alpha(nBodies);
updaters.state_->extractAlpha(alpha);
+
+ Vector totalVelocityIterate;
+ nBodyAssembler.nodalToTransformed(velocityIterates, totalVelocityIterate);
+
for (fixedPointIteration = 0; fixedPointIteration < fixedPointMaxIterations;
++fixedPointIteration) {
// contribution from nonlinearity
globalFriction.updateAlpha(alpha);
- Vector totalVelocityIterate;
- nBodyAssembler.nodalToTransformed(velocityIterates, totalVelocityIterate);
-
solveProblem(contactNetwork, bilinearForm, totalRhs, totalVelocityIterate, totalDirichletNodes, lower, upper);
nBodyAssembler.postprocess(totalVelocityIterate, velocityIterates);
- std::vector<Vector> v_m;
+ std::vector<Vector> v_m; //TODO : wrong, isnt used atm;
updaters.rate_->extractOldVelocity(v_m);
for (size_t i=0; i<v_m.size(); i++) {
@@ -730,6 +739,7 @@ int main(int argc, char *argv[]) { try {
//const auto& stateEnergyNorms = contactNetwork.stateEnergyNorms();
+ timeSteps = parset.get<size_t>("timeSteps.timeSteps");
for (timeStep=1; timeStep<=timeSteps; timeStep++) {
advanceStep(current, contactNetwork);
diff --git a/src/spatial-solving/solverfactory.cc b/src/spatial-solving/solverfactory.cc
index 5f02147e..0650c522 100644
--- a/src/spatial-solving/solverfactory.cc
+++ b/src/spatial-solving/solverfactory.cc
@@ -20,14 +20,10 @@ SolverFactory<Functional, BitVector>::SolverFactory(
J_(Dune::Solvers::wrap_own_share<const Functional>(std::forward<Functional>(J))) {
auto localSolver = Dune::TNNMG::gaussSeidelLocalSolver(LocalSolver());
- //nonlinearSmoother_ = std::make_shared<NonlinearSmoother>(*J_, dummyIterate_, localSolver);
-
- //using MyNonlinearSmoother = Dune::TNNMG::NonlinearGSStep<Functional, decltype(localSolver), BitVector>;
nonlinearSmoother_ = std::make_shared<NonlinearSmoother>(*J_, dummyIterate_, localSolver);
//nonlinearSmoother_ = std::make_shared<NonlinearSmoother>(*J_, dummyIterate_, LocalSolver());
-
auto linearSolver_ptr = Dune::Solvers::wrap_own_share<std::decay_t<LinearSolver>>(std::forward<LinearSolver>(linearSolver));
//tnnmgStep_ = std::make_shared<Step>(*J_, dummyIterate_, nonlinearSmoother_, linearSolver_ptr, 10, DefectProjection(), LineSearchSolver());
diff --git a/src/spatial-solving/solverfactory.hh b/src/spatial-solving/solverfactory.hh
index 5fe7e62f..203bd955 100644
--- a/src/spatial-solving/solverfactory.hh
+++ b/src/spatial-solving/solverfactory.hh
@@ -33,7 +33,7 @@ class SolverFactory {
using Linearization = Linearization<Functional, BitVector>;
using DefectProjection = typename Dune::TNNMG::ObstacleDefectProjection;
- //using Step = TNNMGStep<Functional, BitVector, Linearization, DefectProjection, LineSearchSolver>;
+ //using Step = Dune::TNNMG::TNNMGStep<Functional, BitVector, Linearization, DefectProjection, LineSearchSolver>;
using Step = Dune::TNNMG::TNNMGStep<Functional, BitVector, Linearization, DefectProjection, LineSearchSolver>;
template <class LinearSolver>
--
GitLab