diff --git a/AMDiS/src/nonlin/Newton.h b/AMDiS/src/nonlin/Newton.h new file mode 100644 index 0000000000000000000000000000000000000000..2be9aaeb5d3a118e6aba5b2fed0fddeb4e11fe00 --- /dev/null +++ b/AMDiS/src/nonlin/Newton.h @@ -0,0 +1,161 @@ +// ============================================================================ +// == == +// == 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. + + + +/** \file Newton.h */ + +#ifndef AMDIS_NEWTON_H +#define AMDIS_NEWTON_H + +#include "CreatorInterface.h" +#include "NonLinSolver.h" +#include "OEMSolver.h" + +namespace AMDiS { + + /** + * \ingroup Solver + * + * \Brief + * Implements the newton method for solving a non linear system. Sub class of + * NonLinSolver. + */ + class Newton : public NonLinSolver + { + public: + /// Creator class used in the NonLinSolverMap. + class Creator : public NonLinSolverCreator + { + public: + virtual ~Creator() {} + + /** \brief + * Returns a new Newton object. + */ + NonLinSolver* create() + { + return new Newton(this->name, this->linearSolver, this->nonLinUpdater); + } + }; + + /** \brief + * Calls constructor of base class NonLinSolver + */ + Newton(const std::string& name_, + OEMSolver *linSolver_, + NonLinUpdater *updater) + : NonLinSolver(name_, linSolver_, updater), + b(NULL) + {} + + private: + /** \brief + * realisation of NonLinSolver::init + */ + void init() + {} + + /** \brief + * realisation of NonLinSolver::nlsolve + */ + int nlsolve(SolverMatrix<Matrix<DOFMatrix*> >& mat, + SystemVector& x, SystemVector& rhs) + { + //DOFVector<T> *d = problem->getRHS(); + //DOFVector<T> *x = problem->getSolution();; + + b = new SystemVector(x); + *b = rhs; + + // // copy operators from fh to b + // std::vector<Operator*>::iterator op; + // std::vector<double*>::iterator fac; + // for(op = d->getOperatorsBegin(), + // fac = d->getOperatorFactorBegin(); + // op != d->getOperatorsEnd(); + // ++op, ++fac) + // { + // b->addOperator(*op, *fac); + // } + + double err = 0.0, err_old = -1.0; + int iter, n; + + INFO(this->info,2)("iter. | this->residual | red. | n |\n"); + + for (iter = 1; iter <= this->maxIter; iter++) { + /*--- Assemble DF(x) and F(x) ----------------------------------------------*/ + this->nonLinUpdater->update(/*x,*/true, b); + /*--- Initial guess is zero ------------------------------------------------*/ + rhs.set(0.0); + /*--- solve linear system --------------------------------------------------*/ + n = solveLinearSystem(mat, *b, rhs); + /*--- x = x - d ------------------------------------------------------------*/ + x -= rhs; + + if (this->usedNorm == NO_NORM || this->usedNorm == L2_NORM) + err = L2Norm(&rhs); // sollte hier nicht b genommen werden (s. NewtonS) ? + else + err = H1Norm(&rhs); // sollte hier nicht b genommen werden (s. NewtonS) ? + + + if (iter == 1) this->initial_residual = err; + + if (err_old <= 0) { + INFO(this->info,2)("%5d | %12.5e | -------- | %4d |\n", iter, err, n); + } else { + INFO(this->info,2)("%5d | %12.5e | %8.2e | %4d |\n", + iter, err, err/err_old, n); + } + + if ((this->residual = err) < this->tolerance) { + INFO(this->info,4)("finished successfully\n"); + return iter; + } + err_old = err; + } + + if (this->info < 2) + INFO(this->info,1)("iter. %d, residual: %12.5e\n", iter, err); + INFO(this->info,1)("tolerance %e not reached\n", this->tolerance); + + this->residual = err; + + return iter; + } + + /** \brief + * realisation of NonLinSolver::exit + */ + void exit() + { + if (b != NULL) + delete b; + } + + private: + /** \brief + * internal used data + */ + SystemVector *b; + }; + +} + +#endif // AMDIS_NEWTON_H diff --git a/AMDiS/src/nonlin/NonLinSolver.h b/AMDiS/src/nonlin/NonLinSolver.h index 671a1d1cc03b5b60f9bc73a88e3ed554fc68e14b..16ee8a2875684f4c88af0b757f494029516c3732 100644 --- a/AMDiS/src/nonlin/NonLinSolver.h +++ b/AMDiS/src/nonlin/NonLinSolver.h @@ -181,5 +181,7 @@ namespace AMDiS { } +#include "Newton.h" + #endif // AMDIS_NONLINSOLVER_H