diff --git a/dumux/nonlinear/newtonsolver.hh b/dumux/nonlinear/newtonsolver.hh
index c130bb3f409a5e5f1a3268f57e383e5dbe9614d3..b6a46ae0d2f36381b180b50dbb29472ab7c2f336 100644
--- a/dumux/nonlinear/newtonsolver.hh
+++ b/dumux/nonlinear/newtonsolver.hh
@@ -205,20 +205,6 @@ public:
void setMaxSteps(int maxSteps)
{ maxSteps_ = maxSteps; }
- /*!
- * \brief Run the Newton method to solve a non-linear system.
- * Does time step control when the Newton fails to converge
- */
- [[deprecated("Use attachConvergenceWriter(convWriter) and solve(x, *timeLoop) instead")]]
- void solve(SolutionVector& uCurrentIter, TimeLoop& timeLoop,
- std::shared_ptr<ConvergenceWriter> convWriter)
- {
- if (!convergenceWriter_)
- attachConvergenceWriter(convWriter);
-
- solve(uCurrentIter, timeLoop);
- }
-
/*!
* \brief Run the Newton method to solve a non-linear system.
* Does time step control when the Newton fails to converge
@@ -684,20 +670,6 @@ public:
return oldTimeStep*(1.0 + percent/1.2);
}
- /*!
- * \brief Specifies if the Newton method ought to be chatty.
- */
- [[deprecated("Has been replaced by setVerbosity(int). Will be removed after 3.1 release!")]]
- void setVerbose(bool val)
- { verbosity_ = val; }
-
- /*!
- * \brief Returns true if the Newton method ought to be chatty.
- */
- [[deprecated("Has been replaced by int verbosity(). Will be removed after 3.1 release!")]]
- bool verbose() const
- { return verbosity_ ; }
-
/*!
* \brief Specifies the verbosity level
*/
@@ -992,8 +964,9 @@ private:
void newtonUpdateShiftImpl_(const SolVec &uLastIter,
const SolVec &deltaU)
{
- for (int i = 0; i < int(uLastIter.size()); ++i) {
- typename SolVec::block_type uNewI = uLastIter[i];
+ for (int i = 0; i < int(uLastIter.size()); ++i)
+ {
+ auto uNewI = uLastIter[i];
uNewI -= deltaU[i];
Scalar shiftAtDof = relativeShiftAtDof_(uLastIter[i], uNewI);
@@ -1023,7 +996,6 @@ private:
const SolutionVector &deltaU)
{
Scalar lambda = 1.0;
- SolutionVector tmp(uLastIter);
while (true)
{
@@ -1082,7 +1054,7 @@ private:
//! to this field vector type in Dune ISTL
//! Could be avoided for vectors that already have the right type using SFINAE
//! but it shouldn't impact performance too much
- constexpr auto blockSize = JacobianMatrix::block_type::rows;
+ constexpr auto blockSize = std::decay_t<decltype(b[0])>::dimension;
using BlockType = Dune::FieldVector<Scalar, blockSize>;
Dune::BlockVector<BlockType> xTmp; xTmp.resize(b.size());
Dune::BlockVector<BlockType> bTmp(xTmp);
diff --git a/test/nonlinear/CMakeLists.txt b/test/nonlinear/CMakeLists.txt
index 3640235bdcfbdd1b4fa0aacce2527f8cc8114bf9..7036b25e084f12339a751747c7321c19bf38b612 100644
--- a/test/nonlinear/CMakeLists.txt
+++ b/test/nonlinear/CMakeLists.txt
@@ -1 +1,2 @@
add_subdirectory(findscalarroot)
+add_subdirectory(newton)
diff --git a/test/nonlinear/newton/CMakeLists.txt b/test/nonlinear/newton/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..4187ae5cab66a15cfa3ad925070af87a60174f10
--- /dev/null
+++ b/test/nonlinear/newton/CMakeLists.txt
@@ -0,0 +1,7 @@
+dumux_add_test(SOURCES test_newton.cc
+ LABELS unit nonlinear)
+dumux_add_test(NAME test_newton_linesearch
+ TARGET test_newton
+ COMMAND test_newton
+ CMD_ARGS "-Newton.UseLineSearch" "true"
+ LABELS unit nonlinear)
diff --git a/test/nonlinear/newton/test_newton.cc b/test/nonlinear/newton/test_newton.cc
new file mode 100644
index 0000000000000000000000000000000000000000..5f37e2c3656bbdaa8f0574888023e0823139cd6b
--- /dev/null
+++ b/test/nonlinear/newton/test_newton.cc
@@ -0,0 +1,163 @@
+#include <config.h>
+
+#include <iostream>
+#include <cmath>
+#include <iomanip>
+
+#include <dune/common/exceptions.hh>
+#include <dune/common/float_cmp.hh>
+#include <dune/common/parallel/mpihelper.hh>
+#include <dumux/nonlinear/newtonsolver.hh>
+
+/*
+
+ This test currently solves a scalar non-linear equation using the
+ Dumux::NewtonSolver. The Mock classes expose which dependencies the
+ current implementation has on different other classes it interacts with.
+ Several dependencies seem unnecessary. In particular, the current
+ implementation is basically hard-coded to double indexable residual vectors
+ with several math operators. A good idea would seem to somehow delegate
+ this dependency to something like a linear algebra backend or at least
+ the assembler. The assembler requires a lot of interface functions which
+ are not always needed. The linear solver interdependency is much better (small).
+
+ This test is to ensure that the dependencies do not grow more in the future.
+
+ */
+
+namespace Dumux {
+
+class MockScalarVariables
+{
+public:
+ static constexpr int dimension = 1;
+ MockScalarVariables() : var_(0.0) {}
+ explicit MockScalarVariables(double v) : var_(v) {}
+ MockScalarVariables& operator-=(const MockScalarVariables& other) { var_ -= other.var_; return *this; }
+ double& operator[] (int i) { return var_; }
+ const double& operator[] (int i) const { return var_; }
+private:
+ double var_;
+};
+
+class MockScalarResidual
+{
+public:
+ MockScalarResidual() : res_(0.0) {}
+ explicit MockScalarResidual(double r) : res_(r) {}
+ MockScalarResidual& operator=(double r) { res_[0] = r; return *this; }
+ MockScalarResidual& operator*=(double a) { res_[0] *= a; return *this; }
+ MockScalarResidual& operator+=(const MockScalarResidual& other) { res_[0] += other.res_[0]; return *this; }
+ MockScalarResidual& operator-=(const MockScalarResidual& other) { res_[0] -= other.res_[0]; return *this; }
+ constexpr std::size_t size() const { return 1; }
+ MockScalarVariables& operator[] (int i) { return res_; }
+ const MockScalarVariables& operator[] (int i) const { return res_; }
+ double two_norm2() const { return res_[0]*res_[0]; }
+private:
+ MockScalarVariables res_;
+};
+
+class MockScalarAssembler
+{
+public:
+ using ResidualType = MockScalarResidual;
+ using Scalar = double;
+ using JacobianMatrix = double;
+
+ void setLinearSystem() {}
+
+ bool isStationaryProblem() { return false; }
+
+ ResidualType prevSol() { return ResidualType(0.0); }
+
+ void resetTimeStep(const ResidualType& sol) {}
+
+ void assembleResidual(const ResidualType& sol)
+ {
+ res_[0][0] = sol[0][0]*sol[0][0] - 5.0;
+ }
+
+ void assembleJacobianAndResidual (const ResidualType& sol)
+ {
+ assembleResidual(sol);
+ jac_ = 2.0*sol[0][0];
+ }
+
+ JacobianMatrix& jacobian() { return jac_; }
+
+ ResidualType& residual() { return res_; }
+
+ double residualNorm(const ResidualType& sol)
+ {
+ assembleResidual(sol);
+ return res_[0][0];
+ }
+
+ void updateGridVariables(const ResidualType& sol) {}
+
+private:
+ JacobianMatrix jac_;
+ ResidualType res_;
+};
+
+class MockScalarLinearSolver
+{
+public:
+ void setResidualReduction(double residualReduction) {}
+
+ template<class Vector>
+ bool solve (const double& A, Vector& x, const Vector& b) const
+ {
+ x[0][0] = b[0][0]/A;
+ return true;
+ }
+};
+
+} // end namespace Dumux
+
+int main(int argc, char* argv[]) try
+{
+ using namespace Dumux;
+
+ // maybe initialize MPI
+ Dune::MPIHelper::instance(argc, argv);
+
+ // initialize parameters
+ // TODO this is necessary because there are some global default used in the Newton solver
+ // Do we really need them to be global defaults???
+ Parameters::init(argc, argv);
+
+ // use the Newton solver to find a solution to a scalar equation
+ using Assembler = MockScalarAssembler;
+ using LinearSolver = MockScalarLinearSolver;
+ using Solver = NewtonSolver<Assembler, LinearSolver, DefaultPartialReassembler>;
+
+ auto assembler = std::make_shared<Assembler>();
+ auto linearSolver = std::make_shared<LinearSolver>();
+ auto solver = std::make_shared<Solver>(assembler, linearSolver);
+
+ double initialGuess = 0.1;
+ MockScalarResidual x(initialGuess);
+
+ std::cout << "Solving: x^2 - 5 = 0" << std::endl;
+ solver->solve(x);
+ std::cout << "Solution: " << std::setprecision(15) << x[0][0]
+ << ", exact: " << std::sqrt(5.0)
+ << ", error: " << std::abs(x[0][0]-std::sqrt(5.0))/std::sqrt(5.0)*100 << "%" << std::endl;
+
+ if (Dune::FloatCmp::ne(x[0][0], std::sqrt(5.0), 1e-13))
+ DUNE_THROW(Dune::Exception, "Didn't find correct root: " << std::setprecision(15) << x[0][0] << ", exact: " << std::sqrt(5.0));
+
+ return 0;
+
+}
+catch (const Dune::Exception& e)
+{
+ std::cout << e << std::endl;
+ return 1;
+}
+catch (...)
+{
+ std::cout << "Unknown exception thrown!" << std::endl;
+ return 1;
+}