From bb8ffa361a47f9a4b0b60d984a4cbca51d47888f Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Thu, 31 Jan 2013 12:07:38 +0100
Subject: [PATCH] Be smarter about damping
---
src/one-body-sample.cc | 28 ++++++++++++++++++++++------
src/one-body-sample.parset | 3 ++-
2 files changed, 24 insertions(+), 7 deletions(-)
diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index 77520b29..e25b11c3 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -134,6 +134,13 @@ template <class FunctionMap> void initPython(FunctionMap &functions) {
.toC<typename FunctionMap::Base>(functions);
}
+template <class VectorType>
+double diff_two_norm(VectorType const &v1, VectorType const &v2) {
+ VectorType tmp = v1;
+ tmp -= v2;
+ return tmp.two_norm();
+}
+
int main(int argc, char *argv[]) {
try {
Dune::Timer timer;
@@ -343,11 +350,12 @@ int main(int argc, char *argv[]) {
};
VectorType ud = ud_initial;
+ SingletonVectorType alpha = alpha_initial;
auto const state_fpi_max =
parset.get<size_t>("solver.tnnmg.fixed_point_iterations");
for (size_t run = 1; run <= timesteps; ++run) {
VectorType u;
- SingletonVectorType alpha;
+ double lastCorrection;
stateUpdater->nextTimeStep();
timeSteppingScheme->nextTimeStep();
@@ -393,15 +401,23 @@ int main(int argc, char *argv[]) {
double const fixedPointTolerance =
parset.get<double>("solver.tnnmg.fixed_point_tolerance");
double const damping = parset.get<double>("solver.damping");
+ double const minimalCorrectionReduction =
+ parset.get<double>("solver.minimal_correction_reduction");
+
for (size_t state_fpi = 1; state_fpi <= state_fpi_max; ++state_fpi) {
stateUpdater->solve(ud);
- if (state_fpi == 1)
- stateUpdater->extractState(alpha);
- else {
+ {
SingletonVectorType computed_state;
stateUpdater->extractState(computed_state);
- alpha *= damping;
- alpha.axpy(1.0 - damping, computed_state);
+ double const correction = diff_two_norm(computed_state, alpha);
+ if (state_fpi <= 2 // Let the first two steps pass through unchanged
+ || correction < minimalCorrectionReduction * lastCorrection)
+ alpha = computed_state;
+ else {
+ alpha *= damping;
+ alpha.axpy(1.0 - damping, computed_state);
+ }
+ lastCorrection = correction;
}
solveDisplacementProblem(problem_iterate, alpha);
diff --git a/src/one-body-sample.parset b/src/one-body-sample.parset
index 5b6a7433..aa9c90a7 100644
--- a/src/one-body-sample.parset
+++ b/src/one-body-sample.parset
@@ -23,7 +23,8 @@ width = 5
[solver]
tolerance = 1e-10
-damping = 0.0
+damping = 0.5
+minimal_correction_reduction = 0.5
[solver.tnnmg]
maxiterations = 1000000
--
GitLab