diff --git a/AMDiS/src/nonlin/Newton.h b/AMDiS/src/nonlin/Newton.h
index 7cf430aba217222efec1cfed081b989fb93d7214..7b4f46a1c1876102cc3489d4ad8291ed2ccd62f5 100644
--- a/AMDiS/src/nonlin/Newton.h
+++ b/AMDiS/src/nonlin/Newton.h
@@ -56,8 +56,12 @@ namespace AMDiS {
/// Calls constructor of base class NonLinSolver
Newton(const std::string& name, OEMSolver *linSolver)
: NonLinSolver(name, linSolver),
- b(NULL)
- {}
+ b(NULL),
+ buildCycle(1)
+ {
+
+ Parameters::get(name + "->build cycle", buildCycle);
+ }
private:
/// Realisation of NonLinSolver::init
@@ -80,8 +84,14 @@ namespace AMDiS {
MSG("iter. | this->residual | red. | n |\n");
for (iter = 1; iter <= this->maxIter; iter++) {
- // Assemble DF(x) and F(x)
- prob->buildAfterCoarsen(adaptInfo, 0, true, true);
+ if (iter == 1 || (buildCycle > 0 && (iter-1) % buildCycle == 0)) {
+ this->getLinearSolver()->setMultipleRhs(false);
+ // Assemble DF(x) and F(x)
+ prob->buildAfterCoarsen(adaptInfo, 0, true, true);
+ } else {
+ this->getLinearSolver()->setMultipleRhs(true);
+ prob->buildAfterCoarsen(adaptInfo, 0, false, true);
+ }
// Initial guess is zero
b->set(0.0);
@@ -134,6 +144,12 @@ namespace AMDiS {
private:
/// Internal used data
SystemVector *b;
+
+ /// build matrix every ith iteration,
+ /// 0...build matrix only once,
+ /// i>=1...rebuild matrix in ith solver iteration,
+ /// standard = 1
+ int buildCycle;
};
}
diff --git a/AMDiS/src/nonlin/NewtonArmijo.h b/AMDiS/src/nonlin/NewtonArmijo.h
new file mode 100644
index 0000000000000000000000000000000000000000..451a10825d944f8e8171a6d672565f9d2d53b741
--- /dev/null
+++ b/AMDiS/src/nonlin/NewtonArmijo.h
@@ -0,0 +1,186 @@
+// ============================================================================
+// == ==
+// == AMDiS - Adaptive multidimensional simulations ==
+// == ==
+// == http://www.amdis-fem.org ==
+// == ==
+// ============================================================================
+//
+// Software License for AMDiS
+//
+// Copyright (c) 2010 Dresden University of Technology
+// All rights reserved.
+// Authors: Simon Vey, Thomas Witkowski et al.
+//
+// This file is part of AMDiS
+//
+// See also license.opensource.txt in the distribution.
+
+
+/// literature: some ideas based on "Zur Numerik nichtlinearer Gelichungssysteme (Teil 2)", Werner Vogt, TU Ilmenau, 2004
+
+/** \file NewtonArmijo.h */
+
+#ifndef AMDIS_NEWTON_ARMIJO_H
+#define AMDIS_NEWTON_ARMIJO_H
+
+#include "CreatorInterface.h"
+#include "NonLinSolver.h"
+#include "OEMSolver.h"
+#include "io/VtkWriter.h"
+
+namespace AMDiS {
+
+ /**
+ * \ingroup Solver
+ *
+ * \Brief
+ * Implements the newton method with Armijo-rule for stepsize control
+ * for solving a non linear system. Sub class of NonLinSolver.
+ */
+ class NewtonArmijo : public NonLinSolver
+ {
+ public:
+ /// Creator class used in the NonLinSolverMap.
+ class Creator : public NonLinSolverCreator
+ {
+ public:
+ virtual ~Creator() {}
+
+ /// Returns a new Newton object.
+ NonLinSolver* create()
+ {
+ return new NewtonArmijo(this->name, this->linearSolver);
+ }
+ };
+
+ /// Calls constructor of base class NonLinSolver
+ NewtonArmijo(const std::string& name, OEMSolver *linSolver)
+ : NonLinSolver(name, linSolver),
+ b(NULL),
+ buildCycle(1),
+ delta(1.e-2), // Abstiegsregulator, z.B. 1.e-2, 1.e-4
+ alpha(0.5), // Daempfungsfaktor, z.B. 0.5
+ nLineSearch(5) // maximale Anzahl line-seach Schritte
+ {
+
+ Parameters::get(name + "->build cycle", buildCycle);
+ Parameters::get(name + "->armijo->delta", delta);
+ Parameters::get(name + "->armijo->alpha", alpha);
+ }
+
+ private:
+ /// Realisation of NonLinSolver::init
+ void init() {}
+
+ /// realisation of NonLinSolver::nlsolve
+ int nlsolve(SolverMatrix<Matrix<DOFMatrix*> >& mat,
+ SystemVector& x, SystemVector& rhs,
+ AdaptInfo *adaptInfo,
+ ProblemStat *prob)
+ {
+ FUNCNAME("Newton::nlsolve()");
+
+ if (b == NULL)
+ b = new SystemVector(x);
+
+ double err = 0.0, errOld = -1.0, lambda = 1.0;
+ int iter, n;
+ SystemVector x_test(x);
+
+ MSG("iter. | this->residual | red. | n | lambda |\n");
+
+ for (iter = 1; iter <= this->maxIter; iter++) {
+ if (iter == 1 || (buildCycle > 0 && (iter-1) % buildCycle == 0)) {
+ this->getLinearSolver()->setMultipleRhs(false);
+ // Assemble DF(x) and F(x)
+ prob->buildAfterCoarsen(adaptInfo, 0, true, true);
+ } else {
+ this->getLinearSolver()->setMultipleRhs(true);
+ prob->buildAfterCoarsen(adaptInfo, 0, false, true);
+ }
+
+ // Initial guess is zero
+ b->set(0.0);
+
+ // Solve linear system
+ n = solveLinearSystem(mat, *b, rhs);
+
+ lambda = std::min(1.0, lambda/alpha);
+
+ // x = x + d
+ x_test = x + lambda * (*b);
+
+ for (int k = 0; k < nLineSearch; ++k) {
+ if (this->usedNorm == NO_NORM || this->usedNorm == L2_NORM)
+ err = L2Norm(b);
+ else
+ err = H1Norm(b);
+ // armijo rule
+ if (err <= (1.0 - 2.0 * delta * lambda) * errOld)
+ break;
+
+ lambda *= alpha;
+ x_test = x + lambda * (*b);
+ }
+ x = x_test;
+
+ if (iter == 1)
+ this->initialResidual = err;
+
+ if (errOld <= 0)
+ MSG("%5d | %12.5e | -------- | %4d | %5.2 |\n", iter, err, n, lambda);
+ else
+ MSG("%5d | %12.5e | %8.2e | %4d | %5.2 |\n", iter, err, err/errOld, n,lambda);
+
+ residual = err;
+ if (err < this->tolerance) {
+ MSG("Finished successfully!\n");
+ return iter;
+ }
+ errOld = err;
+ }
+
+ MSG("iter. %d, residual: %12.5e\n", iter, err);
+ MSG("tolerance %e not reached\n", this->tolerance);
+
+ this->residual = err;
+
+ return iter;
+ }
+
+ /// Realisation of NonLinSolver::exit
+ void exit()
+ {
+ if (b != NULL) {
+ delete b;
+ b = NULL;
+ }
+ }
+
+ private:
+ /// Internal used data
+ SystemVector *b;
+
+ /// build matrix every ith iteration,
+ /// 0...build matrix only once,
+ /// i>=1...rebuild matrix in ith solver iteration,
+ /// standard = 1
+ int buildCycle;
+
+ /// Abstiegsregulator, z.B. 1.e-2, 1.e-4
+ /// hinreichender Abstieg bei |f(x_{k+1]})|^2 < (1-2*delta*lambda_k)*|f(x_k})|^2
+ double delta;
+
+ /// Skalierungsfaktor, z.B. 0.5
+ /// Anpassung des Daempfungsfaktors lambda um Faktor alpha:
+ /// lambda_{k+1} = min(lambda_k/alpha, 1)
+ double alpha;
+
+ /// maximale Anzahl line-seach Schritte
+ int nLineSearch;
+ };
+
+}
+
+#endif // AMDIS_NEWTON_ARMIJO_H
diff --git a/AMDiS/src/nonlin/NonLinSolver.h b/AMDiS/src/nonlin/NonLinSolver.h
index d22be7160a9149fffdd36d3e7d69d83303a6d5e2..f865c6f9398f94a9a20f903baef21d97d80f0aa0 100644
--- a/AMDiS/src/nonlin/NonLinSolver.h
+++ b/AMDiS/src/nonlin/NonLinSolver.h
@@ -175,6 +175,7 @@ namespace AMDiS {
}
#include "Newton.h"
+#include "NewtonArmijo.h"
#endif // AMDIS_NONLINSOLVER_H