From 81fed9edb42d9f4c7e9f1be8399542d8824daf77 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Thu, 5 Apr 2018 17:55:45 +0200 Subject: [PATCH 01/19] [test][1p][implicit] Use NewtonSolver --- .../1p/implicit/compressible/test_1p.cc | 33 ++++-------------- .../compressible/test_1p_stationary.cc | 9 ++--- .../1p/implicit/incompressible/test_1pfv.cc | 1 - .../pointsources/test_1pfv_pointsources.cc | 34 ++++--------------- .../test_1pfv_pointsources_timedependent.cc | 34 ++++--------------- .../porousmediumflow/1p/implicit/test_1pfv.cc | 34 ++++--------------- .../1p/implicit/test_1pfv_fracture2d3d.cc | 34 ++++--------------- .../1p/implicit/test_1pfv_network1d3d.cc | 34 ++++--------------- .../1p/implicit/test_1pnifv.cc | 34 ++++--------------- 9 files changed, 52 insertions(+), 195 deletions(-) diff --git a/test/porousmediumflow/1p/implicit/compressible/test_1p.cc b/test/porousmediumflow/1p/implicit/compressible/test_1p.cc index 82c70f15ca..0681e6b119 100644 --- a/test/porousmediumflow/1p/implicit/compressible/test_1p.cc +++ b/test/porousmediumflow/1p/implicit/compressible/test_1p.cc @@ -41,9 +41,8 @@ #include #include -#include +#include #include -#include #include @@ -108,7 +107,6 @@ int main(int argc, char** argv) try using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); auto tEnd = getParam("TimeLoop.TEnd"); auto dt = getParam("TimeLoop.DtInitial"); - auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); // intialize the vtk output module @@ -130,9 +128,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::NewtonController; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // set some check points for the time loop timeLoop->setPeriodicCheckPoint(tEnd/10.0); @@ -143,24 +140,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -176,8 +157,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/compressible/test_1p_stationary.cc b/test/porousmediumflow/1p/implicit/compressible/test_1p_stationary.cc index 9426c67764..61a955a939 100644 --- a/test/porousmediumflow/1p/implicit/compressible/test_1p_stationary.cc +++ b/test/porousmediumflow/1p/implicit/compressible/test_1p_stationary.cc @@ -39,9 +39,8 @@ #include #include -#include +#include #include -#include #include @@ -115,10 +114,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); - using NewtonController = Dumux::NewtonController; - auto newtonController = std::make_shared(); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // linearize & solve Dune::Timer timer; diff --git a/test/porousmediumflow/1p/implicit/incompressible/test_1pfv.cc b/test/porousmediumflow/1p/implicit/incompressible/test_1pfv.cc index 3691b9f81d..5380c2f0fa 100644 --- a/test/porousmediumflow/1p/implicit/incompressible/test_1pfv.cc +++ b/test/porousmediumflow/1p/implicit/incompressible/test_1pfv.cc @@ -34,7 +34,6 @@ #include #include -#include #include #include #include diff --git a/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources.cc b/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources.cc index d15462abb6..92abc5b34e 100644 --- a/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources.cc +++ b/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -104,7 +103,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -132,10 +130,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -143,24 +139,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -175,8 +155,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources_timedependent.cc b/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources_timedependent.cc index 077d34ad20..f0fbf2fe90 100644 --- a/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources_timedependent.cc +++ b/test/porousmediumflow/1p/implicit/pointsources/test_1pfv_pointsources_timedependent.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -104,7 +103,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -132,10 +130,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -146,24 +142,8 @@ int main(int argc, char** argv) try // set the time in the problem for implicit Euler scheme problem->setTime(timeLoop->time() + timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -178,8 +158,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/test_1pfv.cc b/test/porousmediumflow/1p/implicit/test_1pfv.cc index 527af70e84..20b487087d 100644 --- a/test/porousmediumflow/1p/implicit/test_1pfv.cc +++ b/test/porousmediumflow/1p/implicit/test_1pfv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -134,7 +133,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -167,10 +165,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -178,24 +174,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -210,8 +190,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/test_1pfv_fracture2d3d.cc b/test/porousmediumflow/1p/implicit/test_1pfv_fracture2d3d.cc index 957420d572..661a2eacd5 100644 --- a/test/porousmediumflow/1p/implicit/test_1pfv_fracture2d3d.cc +++ b/test/porousmediumflow/1p/implicit/test_1pfv_fracture2d3d.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -128,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -156,10 +154,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -167,24 +163,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -199,8 +179,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/test_1pfv_network1d3d.cc b/test/porousmediumflow/1p/implicit/test_1pfv_network1d3d.cc index 6b61ab1ce5..8002c62815 100644 --- a/test/porousmediumflow/1p/implicit/test_1pfv_network1d3d.cc +++ b/test/porousmediumflow/1p/implicit/test_1pfv_network1d3d.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -128,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -156,10 +154,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -167,24 +163,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // output l2 norm for convergence analysis problem->outputL2Norm(x); @@ -202,8 +182,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1p/implicit/test_1pnifv.cc b/test/porousmediumflow/1p/implicit/test_1pnifv.cc index f9f83626eb..cda24e587c 100644 --- a/test/porousmediumflow/1p/implicit/test_1pnifv.cc +++ b/test/porousmediumflow/1p/implicit/test_1pnifv.cc @@ -42,8 +42,7 @@ #include #include -#include -#include +#include #include #include @@ -130,7 +129,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -162,10 +160,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -173,24 +169,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // compute the new analytical temperature field for the output problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -205,8 +185,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); if (timeLoop->timeStepIndex()==0 || timeLoop->timeStepIndex() % vtkOutputInterval == 0 || timeLoop->willBeFinished()) vtkWriter.write(timeLoop->time()); -- GitLab From e55d8543ada022cf5f88378d1cb69144c1967d68 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 10:31:19 +0200 Subject: [PATCH 02/19] [test][1pnc][implicit] Use NewtonSolver --- .../1pnc/implicit/test_1p2c_fv.cc | 2 +- .../implicit/test_1p2cni_conduction_fv.cc | 33 ++++--------------- .../implicit/test_1p2cni_convection_fv.cc | 33 ++++--------------- 3 files changed, 15 insertions(+), 53 deletions(-) diff --git a/test/porousmediumflow/1pnc/implicit/test_1p2c_fv.cc b/test/porousmediumflow/1pnc/implicit/test_1p2c_fv.cc index 3bd306cfd8..d9bb15ae64 100644 --- a/test/porousmediumflow/1pnc/implicit/test_1p2c_fv.cc +++ b/test/porousmediumflow/1pnc/implicit/test_1p2c_fv.cc @@ -152,7 +152,7 @@ // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt as suggested by the newton solver timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1pnc/implicit/test_1p2cni_conduction_fv.cc b/test/porousmediumflow/1pnc/implicit/test_1p2cni_conduction_fv.cc index 39b9bd5270..d259b57257 100644 --- a/test/porousmediumflow/1pnc/implicit/test_1p2cni_conduction_fv.cc +++ b/test/porousmediumflow/1pnc/implicit/test_1p2cni_conduction_fv.cc @@ -41,9 +41,8 @@ #include #include - #include + #include #include - #include #include @@ -108,7 +107,6 @@ using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); auto tEnd = getParam("TimeLoop.TEnd"); auto dt = getParam("TimeLoop.DtInitial"); - auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); // intialize the vtk output module @@ -134,9 +132,8 @@ auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::NewtonController; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -144,24 +141,8 @@ // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // update the exact time temperature problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -180,8 +161,8 @@ // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/1pnc/implicit/test_1p2cni_convection_fv.cc b/test/porousmediumflow/1pnc/implicit/test_1p2cni_convection_fv.cc index 2d99ef79f9..8bfdfd015d 100644 --- a/test/porousmediumflow/1pnc/implicit/test_1p2cni_convection_fv.cc +++ b/test/porousmediumflow/1pnc/implicit/test_1p2cni_convection_fv.cc @@ -41,9 +41,8 @@ #include #include - #include + #include #include - #include #include @@ -108,7 +107,6 @@ using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); auto tEnd = getParam("TimeLoop.TEnd"); auto dt = getParam("TimeLoop.DtInitial"); - auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); // intialize the vtk output module @@ -134,9 +132,8 @@ auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::NewtonController; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -144,24 +141,8 @@ // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // linearize & solve + nonLinearSolver.solve(x, *timeLoop); // update the exact time temperature problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -180,8 +161,8 @@ // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 45245b211d6b3dde5c9c32a11cafdc403795817b Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 15:47:16 +0200 Subject: [PATCH 03/19] [test][1pncmin][implicit] Use NewtonSolver --- .../1pncmin/implicit/test_1pncminni_fv.cc | 61 +++++++------------ 1 file changed, 21 insertions(+), 40 deletions(-) diff --git a/test/porousmediumflow/1pncmin/implicit/test_1pncminni_fv.cc b/test/porousmediumflow/1pncmin/implicit/test_1pncminni_fv.cc index ad96a94b3d..1c2903a3fe 100644 --- a/test/porousmediumflow/1pncmin/implicit/test_1pncminni_fv.cc +++ b/test/porousmediumflow/1pncmin/implicit/test_1pncminni_fv.cc @@ -30,8 +30,7 @@ #include #include -#include -#include +#include #include #include @@ -134,7 +133,6 @@ int main(int argc, char** argv) try using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); auto tEnd = getParam("TimeLoop.TEnd"); auto dt = getParam("TimeLoop.DtInitial"); - auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); // intialize the vtk output module @@ -164,9 +162,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::NewtonController; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -177,43 +174,27 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } - - // make the new solution the old solution - xOld = x; - gridVariables->advanceTimeStep(); - - // advance to the time loop to the next step - timeLoop->advanceTimeStep(); - - // update the output fields before write - problem->updateVtkOutput(x); + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); + + // make the new solution the old solution + xOld = x; + gridVariables->advanceTimeStep(); + + // advance to the time loop to the next step + timeLoop->advanceTimeStep(); + + // update the output fields before write + problem->updateVtkOutput(x); - // write vtk output - vtkWriter.write(timeLoop->time()); + // write vtk output + vtkWriter.write(timeLoop->time()); - // report statistics of this time step - timeLoop->reportTimeStep(); + // report statistics of this time step + timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 40dfa3414a16fc7822c29a89b41616c93f7d9e55 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 15:52:36 +0200 Subject: [PATCH 04/19] [test][2p][implicit] Use NewtonSolver --- .../implicit/adaptive/test_2p_adaptive_fv.cc | 33 ++++-------------- .../implicit/fracture/test_2p_fracture_fv.cc | 34 ++++--------------- .../2p/implicit/nonisothermal/test_2pni_fv.cc | 34 ++++--------------- 3 files changed, 21 insertions(+), 80 deletions(-) diff --git a/test/porousmediumflow/2p/implicit/adaptive/test_2p_adaptive_fv.cc b/test/porousmediumflow/2p/implicit/adaptive/test_2p_adaptive_fv.cc index dd1c181f56..7ef24730c3 100644 --- a/test/porousmediumflow/2p/implicit/adaptive/test_2p_adaptive_fv.cc +++ b/test/porousmediumflow/2p/implicit/adaptive/test_2p_adaptive_fv.cc @@ -38,8 +38,7 @@ #include #include -#include -#include +#include #include #include @@ -180,7 +179,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -208,10 +206,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -241,23 +237,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -272,8 +253,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/2p/implicit/fracture/test_2p_fracture_fv.cc b/test/porousmediumflow/2p/implicit/fracture/test_2p_fracture_fv.cc index d56429596d..6090dc310d 100644 --- a/test/porousmediumflow/2p/implicit/fracture/test_2p_fracture_fv.cc +++ b/test/porousmediumflow/2p/implicit/fracture/test_2p_fracture_fv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -117,7 +116,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -145,10 +143,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -156,24 +152,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -188,8 +168,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/2p/implicit/nonisothermal/test_2pni_fv.cc b/test/porousmediumflow/2p/implicit/nonisothermal/test_2pni_fv.cc index 1354069095..22b8215f5e 100644 --- a/test/porousmediumflow/2p/implicit/nonisothermal/test_2pni_fv.cc +++ b/test/porousmediumflow/2p/implicit/nonisothermal/test_2pni_fv.cc @@ -38,8 +38,7 @@ #include #include -#include -#include +#include #include @@ -128,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -156,10 +154,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -167,24 +163,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -199,8 +179,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 59a2a1cc280b24b47661aab08063b8c457be61b8 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 15:56:27 +0200 Subject: [PATCH 05/19] [test][2p1cni][implicit] Use NewtonSolver --- .../2p1c/implicit/test_2p1c_fv.cc | 34 +++++-------------- 1 file changed, 8 insertions(+), 26 deletions(-) diff --git a/test/porousmediumflow/2p1c/implicit/test_2p1c_fv.cc b/test/porousmediumflow/2p1c/implicit/test_2p1c_fv.cc index d87098fa99..f6ec0b61d7 100644 --- a/test/porousmediumflow/2p1c/implicit/test_2p1c_fv.cc +++ b/test/porousmediumflow/2p1c/implicit/test_2p1c_fv.cc @@ -41,9 +41,8 @@ #include #include - #include #include - #include + #include #include @@ -108,7 +107,6 @@ using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); auto tEnd = getParam("TimeLoop.TEnd"); auto dt = getParam("TimeLoop.DtInitial"); - auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); // intialize the vtk output module @@ -130,9 +128,9 @@ auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::PriVarSwitchNewtonController; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); + using NewtonSolver = Dumux::PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -140,24 +138,8 @@ // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -172,8 +154,8 @@ // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 78686c9f002ee3fc015a98499e2258cce55308ee Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 15:59:01 +0200 Subject: [PATCH 06/19] [test][2p2c][implicit] Use NewtonSolver --- .../mpnccomparison/test_2p2c_comparison_fv.cc | 35 +++++-------------- .../2p2c/implicit/test_2p2c_fv.cc | 2 +- 2 files changed, 9 insertions(+), 28 deletions(-) diff --git a/test/porousmediumflow/2p2c/implicit/mpnccomparison/test_2p2c_comparison_fv.cc b/test/porousmediumflow/2p2c/implicit/mpnccomparison/test_2p2c_comparison_fv.cc index e8835a9e6d..6d481bae8a 100644 --- a/test/porousmediumflow/2p2c/implicit/mpnccomparison/test_2p2c_comparison_fv.cc +++ b/test/porousmediumflow/2p2c/implicit/mpnccomparison/test_2p2c_comparison_fv.cc @@ -38,8 +38,7 @@ #include #include -#include -#include +#include #include #include @@ -129,7 +128,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -157,10 +155,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -168,24 +165,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -200,8 +181,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/2p2c/implicit/test_2p2c_fv.cc b/test/porousmediumflow/2p2c/implicit/test_2p2c_fv.cc index c7fc2a9a98..1b86e04ba7 100644 --- a/test/porousmediumflow/2p2c/implicit/test_2p2c_fv.cc +++ b/test/porousmediumflow/2p2c/implicit/test_2p2c_fv.cc @@ -153,7 +153,7 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt as suggested by the newton solver timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write vtk output -- GitLab From 80be51d9d659e8560303454f06b77fccf7b08552 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:02:37 +0200 Subject: [PATCH 07/19] [test][2pnc][implicit] Use NewtonSolver --- .../2pnc/implicit/test_2pnc_fv.cc | 35 +++++-------------- .../2pnc/implicit/test_cc2pnc_diffusion.cc | 35 +++++-------------- 2 files changed, 16 insertions(+), 54 deletions(-) diff --git a/test/porousmediumflow/2pnc/implicit/test_2pnc_fv.cc b/test/porousmediumflow/2pnc/implicit/test_2pnc_fv.cc index 9b75b3116f..21750d4981 100644 --- a/test/porousmediumflow/2pnc/implicit/test_2pnc_fv.cc +++ b/test/porousmediumflow/2pnc/implicit/test_2pnc_fv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -128,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -157,10 +155,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -168,24 +165,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -203,8 +184,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/2pnc/implicit/test_cc2pnc_diffusion.cc b/test/porousmediumflow/2pnc/implicit/test_cc2pnc_diffusion.cc index 88ee0bfbbe..31870b972c 100644 --- a/test/porousmediumflow/2pnc/implicit/test_cc2pnc_diffusion.cc +++ b/test/porousmediumflow/2pnc/implicit/test_cc2pnc_diffusion.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -126,7 +125,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -154,10 +152,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -165,24 +162,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -197,8 +178,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From a0883b83ba91db16267faad4f70dd380bbe1c0aa Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:04:03 +0200 Subject: [PATCH 08/19] [test][2pncmin][implicit] Use NewtonSolver --- .../2pncmin/implicit/test_2pncmin_fv.cc | 35 +++++-------------- 1 file changed, 8 insertions(+), 27 deletions(-) diff --git a/test/porousmediumflow/2pncmin/implicit/test_2pncmin_fv.cc b/test/porousmediumflow/2pncmin/implicit/test_2pncmin_fv.cc index f40e7f9e16..1ccc7b74a2 100644 --- a/test/porousmediumflow/2pncmin/implicit/test_2pncmin_fv.cc +++ b/test/porousmediumflow/2pncmin/implicit/test_2pncmin_fv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -126,7 +125,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -158,10 +156,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -173,24 +170,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -208,8 +189,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From e7c537232b8e7307043083fa943c5821d92b11a0 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:07:42 +0200 Subject: [PATCH 09/19] [test][3p][implicit] Use NewtonSolver --- .../3p/implicit/test_3p_fv.cc | 34 ++++--------------- .../3p/implicit/test_3pni_fv_conduction.cc | 34 ++++--------------- .../3p/implicit/test_3pni_fv_convection.cc | 34 ++++--------------- 3 files changed, 21 insertions(+), 81 deletions(-) diff --git a/test/porousmediumflow/3p/implicit/test_3p_fv.cc b/test/porousmediumflow/3p/implicit/test_3p_fv.cc index 8f7a7eda6a..6a17a142b3 100644 --- a/test/porousmediumflow/3p/implicit/test_3p_fv.cc +++ b/test/porousmediumflow/3p/implicit/test_3p_fv.cc @@ -39,8 +39,7 @@ #include #include -#include -#include +#include #include #include @@ -133,7 +132,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -161,10 +159,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -176,24 +172,8 @@ int main(int argc, char** argv) try // the boundary conditions for the implicit Euler scheme problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -205,8 +185,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write vtk output vtkWriter.write(timeLoop->time()); diff --git a/test/porousmediumflow/3p/implicit/test_3pni_fv_conduction.cc b/test/porousmediumflow/3p/implicit/test_3pni_fv_conduction.cc index d57c83212c..7ea429bfed 100644 --- a/test/porousmediumflow/3p/implicit/test_3pni_fv_conduction.cc +++ b/test/porousmediumflow/3p/implicit/test_3pni_fv_conduction.cc @@ -39,8 +39,7 @@ #include #include -#include -#include +#include #include #include @@ -133,7 +132,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -164,10 +162,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -175,24 +171,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // update the exact time temperature problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -207,8 +187,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); if (timeLoop->timeStepIndex()==0 || timeLoop->timeStepIndex() % vtkOutputInterval == 0 || timeLoop->willBeFinished()) vtkWriter.write(timeLoop->time()); diff --git a/test/porousmediumflow/3p/implicit/test_3pni_fv_convection.cc b/test/porousmediumflow/3p/implicit/test_3pni_fv_convection.cc index b608af5e75..dcaf2bcac8 100644 --- a/test/porousmediumflow/3p/implicit/test_3pni_fv_convection.cc +++ b/test/porousmediumflow/3p/implicit/test_3pni_fv_convection.cc @@ -39,8 +39,7 @@ #include #include -#include -#include +#include #include #include @@ -133,7 +132,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -164,10 +162,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -175,24 +171,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // compute the new analytical temperature field for the output problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -207,8 +187,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write vtk output if(timeLoop->timeStepIndex()==0 || timeLoop->timeStepIndex() % vtkOutputInterval == 0 || timeLoop->willBeFinished()) -- GitLab From 021e3c26999e9527cc34f628d330a211d363eee7 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:09:14 +0200 Subject: [PATCH 10/19] [test][3p3c][implicit] Use NewtonSolver --- .../3p3c/implicit/test_3p3c_fv.cc | 35 +++++-------------- 1 file changed, 8 insertions(+), 27 deletions(-) diff --git a/test/porousmediumflow/3p3c/implicit/test_3p3c_fv.cc b/test/porousmediumflow/3p3c/implicit/test_3p3c_fv.cc index a9d8a49efb..d656cad00b 100644 --- a/test/porousmediumflow/3p3c/implicit/test_3p3c_fv.cc +++ b/test/porousmediumflow/3p3c/implicit/test_3p3c_fv.cc @@ -37,8 +37,7 @@ #include #include -#include -#include +#include #include #include @@ -133,7 +132,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -162,10 +160,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -174,24 +171,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -208,8 +189,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } timeLoop->finalize(leafGridView.comm()); -- GitLab From dbfa6b00491e6b782faaf8906e871f70463637a8 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:11:08 +0200 Subject: [PATCH 11/19] [test][3pwateroil][implicit] Use NewtonSolver (fix docu) --- test/porousmediumflow/3pwateroil/implicit/test_box3pwateroil.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/porousmediumflow/3pwateroil/implicit/test_box3pwateroil.cc b/test/porousmediumflow/3pwateroil/implicit/test_box3pwateroil.cc index fd42d37384..bbf72d2633 100644 --- a/test/porousmediumflow/3pwateroil/implicit/test_box3pwateroil.cc +++ b/test/porousmediumflow/3pwateroil/implicit/test_box3pwateroil.cc @@ -180,7 +180,7 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt as suggested by the newton solver timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write vtk output -- GitLab From 870f528eb57282c696609a8899ec300174a93f80 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 16:11:39 +0200 Subject: [PATCH 12/19] [test][co2][implicit] Use NewtonSolver --- .../co2/implicit/test_co2_fv.cc | 35 +++++-------------- 1 file changed, 8 insertions(+), 27 deletions(-) diff --git a/test/porousmediumflow/co2/implicit/test_co2_fv.cc b/test/porousmediumflow/co2/implicit/test_co2_fv.cc index 88162d0d80..c9e2297ebe 100644 --- a/test/porousmediumflow/co2/implicit/test_co2_fv.cc +++ b/test/porousmediumflow/co2/implicit/test_co2_fv.cc @@ -38,8 +38,7 @@ #include #include -#include -#include +#include #include #include @@ -103,7 +102,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -132,10 +130,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -143,24 +140,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -172,8 +153,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write vtk output vtkWriter.write(timeLoop->time()); -- GitLab From 523ce5cfb6eb7712a00195162f5fcb09a0697444 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 17:05:46 +0200 Subject: [PATCH 13/19] [test][richards][implicit] Use NewtonSolver --- .../implicit/test_ccrichardsanalytical.cc | 35 ++++-------------- .../richards/implicit/test_richardslens_fv.cc | 2 +- .../implicit/test_richardsniconduction_fv.cc | 36 +++++-------------- .../implicit/test_richardsniconvection_fv.cc | 36 +++++-------------- .../implicit/test_richardsnievaporation_fv.cc | 2 +- 5 files changed, 25 insertions(+), 86 deletions(-) diff --git a/test/porousmediumflow/richards/implicit/test_ccrichardsanalytical.cc b/test/porousmediumflow/richards/implicit/test_ccrichardsanalytical.cc index e5c0bd42a0..8b9be9df43 100644 --- a/test/porousmediumflow/richards/implicit/test_ccrichardsanalytical.cc +++ b/test/porousmediumflow/richards/implicit/test_ccrichardsanalytical.cc @@ -41,9 +41,7 @@ #include #include -#include -#include -#include +#include #include @@ -129,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -157,10 +154,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::RichardsNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::RichardsNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -169,24 +164,8 @@ int main(int argc, char** argv) try assembler->setPreviousSolution(xOld); problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -203,8 +182,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/richards/implicit/test_richardslens_fv.cc b/test/porousmediumflow/richards/implicit/test_richardslens_fv.cc index 9d260e610e..1550ec49da 100644 --- a/test/porousmediumflow/richards/implicit/test_richardslens_fv.cc +++ b/test/porousmediumflow/richards/implicit/test_richardslens_fv.cc @@ -179,7 +179,7 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt as suggested by the newton solver timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/richards/implicit/test_richardsniconduction_fv.cc b/test/porousmediumflow/richards/implicit/test_richardsniconduction_fv.cc index 49a276444d..82f3c629c1 100644 --- a/test/porousmediumflow/richards/implicit/test_richardsniconduction_fv.cc +++ b/test/porousmediumflow/richards/implicit/test_richardsniconduction_fv.cc @@ -41,9 +41,7 @@ #include #include -#include -#include -#include +#include #include @@ -129,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -158,10 +155,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::RichardsNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::RichardsNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -169,24 +164,9 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); + // compute the new analytical temperature field for the output problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -203,8 +183,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/richards/implicit/test_richardsniconvection_fv.cc b/test/porousmediumflow/richards/implicit/test_richardsniconvection_fv.cc index 4b549054c9..ba18b9ce98 100644 --- a/test/porousmediumflow/richards/implicit/test_richardsniconvection_fv.cc +++ b/test/porousmediumflow/richards/implicit/test_richardsniconvection_fv.cc @@ -41,9 +41,7 @@ #include #include -#include -#include -#include +#include #include @@ -129,7 +127,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -158,10 +155,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::RichardsNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::RichardsNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -169,24 +164,9 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); + // compute the new analytical temperature field for the output problem->updateExactTemperature(x, timeLoop->time()+timeLoop->timeStepSize()); @@ -203,8 +183,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/richards/implicit/test_richardsnievaporation_fv.cc b/test/porousmediumflow/richards/implicit/test_richardsnievaporation_fv.cc index bf53af86c7..3464a371f3 100644 --- a/test/porousmediumflow/richards/implicit/test_richardsnievaporation_fv.cc +++ b/test/porousmediumflow/richards/implicit/test_richardsnievaporation_fv.cc @@ -147,7 +147,7 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt as suggested by the newton solver timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From e014e665c965225de6ecb5d96cfc3a56c48cc3e7 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Fri, 6 Apr 2018 17:22:47 +0200 Subject: [PATCH 14/19] [test][richardsnc][implicit] Use NewtonSolver --- dumux/porousmediumflow/richardsnc/model.hh | 4 +++ .../richardsnc/volumevariables.hh | 2 ++ .../richardsnc/implicit/test_richardsnc_fv.cc | 35 ++++--------------- 3 files changed, 13 insertions(+), 28 deletions(-) diff --git a/dumux/porousmediumflow/richardsnc/model.hh b/dumux/porousmediumflow/richardsnc/model.hh index 5dad8ad9cf..659a42feb0 100644 --- a/dumux/porousmediumflow/richardsnc/model.hh +++ b/dumux/porousmediumflow/richardsnc/model.hh @@ -148,6 +148,10 @@ SET_INT_PROP(RichardsNC, ReplaceCompEqIdx, 0); //! define the VolumeVariables SET_TYPE_PROP(RichardsNC, VolumeVariables, RichardsNCVolumeVariables); +//! The default richardsnc model computes no diffusion in the air phase +//! Turning this on leads to the extended Richards equation (see e.g. Vanderborght et al. 2017) +SET_BOOL_PROP(RichardsNC, EnableWaterDiffusionInAir, false); + /*! *\brief The fluid system used by the model. * diff --git a/dumux/porousmediumflow/richardsnc/volumevariables.hh b/dumux/porousmediumflow/richardsnc/volumevariables.hh index 668e212a7a..cb12812843 100644 --- a/dumux/porousmediumflow/richardsnc/volumevariables.hh +++ b/dumux/porousmediumflow/richardsnc/volumevariables.hh @@ -319,6 +319,8 @@ class RichardsNCVolumeVariables : public RichardsBaseVolumeVariables using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem); using Indices = typename GET_PROP_TYPE(TypeTag, Indices); + static_assert(!GET_PROP_VALUE(TypeTag, EnableWaterDiffusionInAir), "Water diffusion in air is not implement for RichardsNC"); + enum { wPhaseIdx = Indices::wPhaseIdx, diff --git a/test/porousmediumflow/richardsnc/implicit/test_richardsnc_fv.cc b/test/porousmediumflow/richardsnc/implicit/test_richardsnc_fv.cc index 8fcd51145e..22decb0f04 100644 --- a/test/porousmediumflow/richardsnc/implicit/test_richardsnc_fv.cc +++ b/test/porousmediumflow/richardsnc/implicit/test_richardsnc_fv.cc @@ -41,9 +41,7 @@ #include #include -#include -#include -#include +#include #include @@ -130,7 +128,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -158,10 +155,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(); // the non-linear solver - using NewtonController = Dumux::RichardsNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::RichardsNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -169,24 +164,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -202,8 +181,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From faff379fed7ec642ba94eb6330bccc5ebcae9880 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Mon, 9 Apr 2018 06:38:22 +0200 Subject: [PATCH 15/19] [test][tracer][implicit] Use NewtonSolver --- .../tracer/1ptracer/test_cctracer.cc | 3 +- .../tracer/constvel/test_tracer.cc | 3 +- .../multicomp/test_tracer_maxwellstefan.cc | 34 ++++--------------- 3 files changed, 9 insertions(+), 31 deletions(-) diff --git a/test/porousmediumflow/tracer/1ptracer/test_cctracer.cc b/test/porousmediumflow/tracer/1ptracer/test_cctracer.cc index 9187fa0257..010fa97380 100644 --- a/test/porousmediumflow/tracer/1ptracer/test_cctracer.cc +++ b/test/porousmediumflow/tracer/1ptracer/test_cctracer.cc @@ -39,7 +39,6 @@ #include #include -#include #include #include @@ -286,7 +285,7 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller + // set new dt timeLoop->setTimeStepSize(dt); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/tracer/constvel/test_tracer.cc b/test/porousmediumflow/tracer/constvel/test_tracer.cc index f75e2b5a3a..b710725758 100644 --- a/test/porousmediumflow/tracer/constvel/test_tracer.cc +++ b/test/porousmediumflow/tracer/constvel/test_tracer.cc @@ -37,7 +37,6 @@ #include #include -#include #include #include @@ -183,7 +182,7 @@ int main(int argc, char** argv) try if (timeLoop->isCheckPoint() || timeLoop->finished()) vtkWriter.write(timeLoop->time()); - // set new dt as suggested by newton controller + // set new dt timeLoop->setTimeStepSize(dt); } diff --git a/test/porousmediumflow/tracer/multicomp/test_tracer_maxwellstefan.cc b/test/porousmediumflow/tracer/multicomp/test_tracer_maxwellstefan.cc index 07d3eb609d..a03b1a28a8 100644 --- a/test/porousmediumflow/tracer/multicomp/test_tracer_maxwellstefan.cc +++ b/test/porousmediumflow/tracer/multicomp/test_tracer_maxwellstefan.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -135,7 +134,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -163,10 +161,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -174,24 +170,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -207,8 +187,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 006b17733d74e4993c0e410a463483d03dd81646 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Mon, 9 Apr 2018 07:00:32 +0200 Subject: [PATCH 16/19] [test][mpnc][implicit] Use NewtonSolver --- .../nonequilibrium/newtonsolver.hh | 58 +++++++++++++++++++ .../2p2ccomparison/test_mpnc_comparison_fv.cc | 34 +++-------- .../mpnc/implicit/test_boxmpnckinetic.cc | 34 +++-------- .../implicit/test_boxmpncthermalnonequil.cc | 34 +++-------- .../mpnc/implicit/test_mpnc_obstacle_fv.cc | 34 +++-------- 5 files changed, 86 insertions(+), 108 deletions(-) create mode 100644 dumux/porousmediumflow/nonequilibrium/newtonsolver.hh diff --git a/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh b/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh new file mode 100644 index 0000000000..f12622d634 --- /dev/null +++ b/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh @@ -0,0 +1,58 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 2 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + *****************************************************************************/ +/*! + * \file + * \ingroup PorousmediumNonEquilibriumModel + * \brief A MpNc specific newton solver. + * This controller calls the velocity averaging in the model after each iteration. + */ +#ifndef DUMUX_NONEQUILIBRIUM_NEWTON_SOLVER_HH +#define DUMUX_NONEQUILIBRIUM_NEWTON_SOLVER_HH + +#include + +namespace Dumux { +/*! + * \ingroup PorousmediumNonEquilibriumModel + * \brief A nonequilibrium specific newton solver. + * This solver calls the velocity averaging in the problem after each iteration. + */ +template +class NonEquilibriumNewtonSolver : public NewtonSolver +{ + using ParentType = NewtonSolver; + using SolutionVector = typename Assembler::ResidualType; + +public: + using ParentType::ParentType; + + void newtonEndStep(SolutionVector &uCurrentIter, + const SolutionVector &uLastIter) final + { + ParentType::newtonEndStep(uCurrentIter, uLastIter); + + auto& gridVariables = this->assembler().gridVariables(); + // Averages the face velocities of a vertex. Implemented in the model. + // The velocities are stored in the model. + gridVariables.calcVelocityAverage(uCurrentIter); + } +}; + +} // end namespace Dumux +#endif // DUMUX_NONEQUILIBRIUM_NEWTON_SOLVER_HH diff --git a/test/porousmediumflow/mpnc/implicit/2p2ccomparison/test_mpnc_comparison_fv.cc b/test/porousmediumflow/mpnc/implicit/2p2ccomparison/test_mpnc_comparison_fv.cc index 6b64b18c15..4b01fdd4f9 100644 --- a/test/porousmediumflow/mpnc/implicit/2p2ccomparison/test_mpnc_comparison_fv.cc +++ b/test/porousmediumflow/mpnc/implicit/2p2ccomparison/test_mpnc_comparison_fv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -135,7 +134,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -163,10 +161,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -174,24 +170,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -206,8 +186,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/mpnc/implicit/test_boxmpnckinetic.cc b/test/porousmediumflow/mpnc/implicit/test_boxmpnckinetic.cc index 8fd4d91754..2f6e9f2ef0 100644 --- a/test/porousmediumflow/mpnc/implicit/test_boxmpnckinetic.cc +++ b/test/porousmediumflow/mpnc/implicit/test_boxmpnckinetic.cc @@ -40,8 +40,7 @@ #include #include -#include -#include +#include #include #include @@ -134,7 +133,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -162,10 +160,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = NonEquilibriumNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NonEquilibriumNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -173,24 +169,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -205,8 +185,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/mpnc/implicit/test_boxmpncthermalnonequil.cc b/test/porousmediumflow/mpnc/implicit/test_boxmpncthermalnonequil.cc index 8974eecafa..7501c698c7 100644 --- a/test/porousmediumflow/mpnc/implicit/test_boxmpncthermalnonequil.cc +++ b/test/porousmediumflow/mpnc/implicit/test_boxmpncthermalnonequil.cc @@ -40,8 +40,7 @@ #include #include -#include -#include +#include #include #include @@ -133,7 +132,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -161,10 +159,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = NonEquilibriumNewtonController; - using NewtonMethod = Dumux::NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NonEquilibriumNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -172,24 +168,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -205,8 +185,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); diff --git a/test/porousmediumflow/mpnc/implicit/test_mpnc_obstacle_fv.cc b/test/porousmediumflow/mpnc/implicit/test_mpnc_obstacle_fv.cc index fc999af2a1..86bdae2d55 100644 --- a/test/porousmediumflow/mpnc/implicit/test_mpnc_obstacle_fv.cc +++ b/test/porousmediumflow/mpnc/implicit/test_mpnc_obstacle_fv.cc @@ -41,8 +41,7 @@ #include #include -#include -#include +#include #include #include @@ -135,7 +134,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -163,10 +161,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -174,24 +170,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -206,8 +186,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); } while (!timeLoop->finished()); -- GitLab From 76995aedd9e54c44844d164329994399c6043cc4 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Mon, 9 Apr 2018 07:13:16 +0200 Subject: [PATCH 17/19] [tutorial] Use NewtonSolver --- tutorial/ex1/exercise1.cc | 34 ++++-------------- tutorial/ex1/exercise1_2p.cc | 34 ++++-------------- tutorial/ex1/exercise1_2p2c.cc | 35 +++++-------------- tutorial/ex2/exercise2.cc | 35 +++++-------------- tutorial/ex3/exercise3.cc | 34 ++++-------------- tutorial/solution/ex1/exercise1_2p2c.cc | 35 +++++-------------- .../solution/ex1/exercise1_2pni_solution.cc | 34 ++++-------------- tutorial/solution/ex2/exercise2_solution.cc | 35 +++++-------------- 8 files changed, 60 insertions(+), 216 deletions(-) diff --git a/tutorial/ex1/exercise1.cc b/tutorial/ex1/exercise1.cc index b7655dac87..d6bb2599a0 100644 --- a/tutorial/ex1/exercise1.cc +++ b/tutorial/ex1/exercise1.cc @@ -36,8 +36,7 @@ #include #include -#include -#include +#include #include #include @@ -103,7 +102,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -125,10 +123,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -140,24 +136,8 @@ int main(int argc, char** argv) try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -169,8 +149,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // output to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/ex1/exercise1_2p.cc b/tutorial/ex1/exercise1_2p.cc index 9cba2e73b2..7ae3bca427 100644 --- a/tutorial/ex1/exercise1_2p.cc +++ b/tutorial/ex1/exercise1_2p.cc @@ -36,8 +36,7 @@ #include #include -#include -#include +#include #include #include @@ -106,7 +105,6 @@ int main(int argc, char** argv) try // of type TYPE given in the group GROUPNAME from the input file using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -128,10 +126,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -143,24 +139,8 @@ int main(int argc, char** argv) try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -172,8 +152,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // output to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/ex1/exercise1_2p2c.cc b/tutorial/ex1/exercise1_2p2c.cc index a0b5ebb035..6f7e515c04 100644 --- a/tutorial/ex1/exercise1_2p2c.cc +++ b/tutorial/ex1/exercise1_2p2c.cc @@ -36,8 +36,7 @@ #include #include -#include -#include +#include #include #include @@ -106,7 +105,6 @@ int main(int argc, char** argv) try // of type TYPE given in the group GROUPNAME from the input file using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -128,10 +126,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); + using NewtonSolver = Dumux::PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -143,24 +140,8 @@ int main(int argc, char** argv) try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -172,8 +153,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // output to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/ex2/exercise2.cc b/tutorial/ex2/exercise2.cc index dba708fbc1..d031e8753b 100644 --- a/tutorial/ex2/exercise2.cc +++ b/tutorial/ex2/exercise2.cc @@ -37,8 +37,7 @@ #include #include -#include -#include +#include #include #include @@ -101,7 +100,6 @@ int main(int argc, char** argv)try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -128,10 +126,9 @@ int main(int argc, char** argv)try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); + using NewtonSolver = Dumux::PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -142,24 +139,8 @@ int main(int argc, char** argv)try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -171,8 +152,8 @@ int main(int argc, char** argv)try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); vtkWriter.write(timeLoop->time()); diff --git a/tutorial/ex3/exercise3.cc b/tutorial/ex3/exercise3.cc index 7111c6a3e0..4669ea7f18 100644 --- a/tutorial/ex3/exercise3.cc +++ b/tutorial/ex3/exercise3.cc @@ -42,8 +42,7 @@ #include #include -#include -#include +#include #include #include @@ -140,7 +139,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -167,10 +165,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -178,24 +174,8 @@ int main(int argc, char** argv) try // set previous solution for storage evaluations assembler->setPreviousSolution(xOld); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -207,8 +187,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // write solution and grid to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/solution/ex1/exercise1_2p2c.cc b/tutorial/solution/ex1/exercise1_2p2c.cc index 313486242d..779a668a61 100644 --- a/tutorial/solution/ex1/exercise1_2p2c.cc +++ b/tutorial/solution/ex1/exercise1_2p2c.cc @@ -36,8 +36,7 @@ #include #include -#include -#include +#include #include #include @@ -103,7 +102,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -125,10 +123,9 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); + using NewtonSolver = Dumux::PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -140,24 +137,8 @@ int main(int argc, char** argv) try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -169,8 +150,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // output to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/solution/ex1/exercise1_2pni_solution.cc b/tutorial/solution/ex1/exercise1_2pni_solution.cc index bd2b3e8357..877fe41e39 100644 --- a/tutorial/solution/ex1/exercise1_2pni_solution.cc +++ b/tutorial/solution/ex1/exercise1_2pni_solution.cc @@ -36,8 +36,7 @@ #include #include -#include -#include +#include #include #include @@ -103,7 +102,6 @@ int main(int argc, char** argv) try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -125,10 +123,8 @@ int main(int argc, char** argv) try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = Dumux::NewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using NewtonSolver = Dumux::NewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); @@ -140,24 +136,8 @@ int main(int argc, char** argv) try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -169,8 +149,8 @@ int main(int argc, char** argv) try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); // output to vtk vtkWriter.write(timeLoop->time()); diff --git a/tutorial/solution/ex2/exercise2_solution.cc b/tutorial/solution/ex2/exercise2_solution.cc index dba708fbc1..d031e8753b 100644 --- a/tutorial/solution/ex2/exercise2_solution.cc +++ b/tutorial/solution/ex2/exercise2_solution.cc @@ -37,8 +37,7 @@ #include #include -#include -#include +#include #include #include @@ -101,7 +100,6 @@ int main(int argc, char** argv)try // get some time loop parameters using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); const auto tEnd = getParam("TimeLoop.TEnd"); - const auto maxDivisions = getParam("TimeLoop.MaxTimeStepDivisions"); const auto maxDt = getParam("TimeLoop.MaxTimeStepSize"); auto dt = getParam("TimeLoop.DtInitial"); @@ -128,10 +126,9 @@ int main(int argc, char** argv)try auto linearSolver = std::make_shared(leafGridView, fvGridGeometry->dofMapper()); // the non-linear solver - using NewtonController = PriVarSwitchNewtonController; - using NewtonMethod = NewtonMethod; - auto newtonController = std::make_shared(timeLoop); - NewtonMethod nonLinearSolver(newtonController, assembler, linearSolver); + using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); + using NewtonSolver = Dumux::PriVarSwitchNewtonSolver; + NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop timeLoop->start(); do @@ -142,24 +139,8 @@ int main(int argc, char** argv)try //set time in problem (is used in time-dependent Neumann boundary condition) problem->setTime(timeLoop->time()+timeLoop->timeStepSize()); - // try solving the non-linear system - for (int i = 0; i < maxDivisions; ++i) - { - // linearize & solve - auto converged = nonLinearSolver.solve(x); - - if (converged) - break; - - if (!converged && i == maxDivisions-1) - DUNE_THROW(Dune::MathError, - "Newton solver didn't converge after " - << maxDivisions - << " time-step divisions. dt=" - << timeLoop->timeStepSize() - << ".\nThe solutions of the current and the previous time steps " - << "have been saved to restart files."); - } + // solve the non-linear system with time step control + nonLinearSolver.solve(x, *timeLoop); // make the new solution the old solution xOld = x; @@ -171,8 +152,8 @@ int main(int argc, char** argv)try // report statistics of this time step timeLoop->reportTimeStep(); - // set new dt as suggested by newton controller - timeLoop->setTimeStepSize(newtonController->suggestTimeStepSize(timeLoop->timeStepSize())); + // set new dt as suggested by the newton solver + timeLoop->setTimeStepSize(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize())); vtkWriter.write(timeLoop->time()); -- GitLab From d99d5076b415a1d2bd517d7330ba15a11fdec9ea Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Mon, 9 Apr 2018 07:26:20 +0200 Subject: [PATCH 18/19] Delete obsolete NewtonController and NewtonMethod headers --- dumux/nonlinear/CMakeLists.txt | 3 +- dumux/nonlinear/newtoncontroller.hh | 651 ------------------ dumux/nonlinear/newtonmethod.hh | 222 ------ dumux/nonlinear/staggerednewtoncontroller.hh | 299 -------- dumux/porousmediumflow/2pnc/CMakeLists.txt | 1 - .../privarswitchnewtoncontroller.hh | 220 ------ .../nonequilibrium/CMakeLists.txt | 2 +- .../nonequilibrium/newtoncontroller.hh | 66 -- .../porousmediumflow/richards/CMakeLists.txt | 3 +- .../richards/newtoncontroller.hh | 122 ---- dumux/porousmediumflow/richardsnc/model.hh | 1 - 11 files changed, 4 insertions(+), 1586 deletions(-) delete mode 100644 dumux/nonlinear/newtoncontroller.hh delete mode 100644 dumux/nonlinear/newtonmethod.hh delete mode 100644 dumux/nonlinear/staggerednewtoncontroller.hh delete mode 100644 dumux/porousmediumflow/compositional/privarswitchnewtoncontroller.hh delete mode 100644 dumux/porousmediumflow/nonequilibrium/newtoncontroller.hh delete mode 100644 dumux/porousmediumflow/richards/newtoncontroller.hh diff --git a/dumux/nonlinear/CMakeLists.txt b/dumux/nonlinear/CMakeLists.txt index 25e858ce85..0bc9109ce1 100644 --- a/dumux/nonlinear/CMakeLists.txt +++ b/dumux/nonlinear/CMakeLists.txt @@ -1,7 +1,6 @@ #install headers install(FILES -newtoncontroller.hh newtonconvergencewriter.hh -newtonmethod.hh +newtonsolver.hh DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/nonlinear) diff --git a/dumux/nonlinear/newtoncontroller.hh b/dumux/nonlinear/newtoncontroller.hh deleted file mode 100644 index a64002a69d..0000000000 --- a/dumux/nonlinear/newtoncontroller.hh +++ /dev/null @@ -1,651 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/**************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup Nonlinear - * \brief Reference implementation of a controller class for the Newton solver. - * - * Usually this controller should be sufficient. - */ -#ifndef DUMUX_NEWTON_CONTROLLER_HH -#define DUMUX_NEWTON_CONTROLLER_HH - -#include - -#include -#include -#include -#include - -#include -#include -#include - -#warning "This file is deprecated. Use NewtonSolver instead." - -namespace Dumux { - -/*! - * \ingroup Nonlinear - * \brief An implementation of a Newton controller - * \tparam Scalar the scalar type - * \tparam Comm the communication object used to communicate with all processes - * \note If you want to specialize only some methods but are happy with the - * defaults of the reference controller, derive your controller from - * this class and simply overload the required methods. - */ -template > -class DUNE_DEPRECATED_MSG("Use NewtonSolver instead.") NewtonController -{ - -public: - //! the communication type used to communicate with all processes - using Communication = Comm; - - /*! - * \brief Constructor for stationary problems - */ - NewtonController(const Communication& comm = Dune::MPIHelper::getCollectiveCommunication(), - const std::string& paramGroup = "") - : comm_(comm) - , endIterMsgStream_(std::ostringstream::out) - , paramGroup_(paramGroup) - { - initParams_(paramGroup); - } - - /*! - * \brief Constructor for instationary problems - */ - NewtonController(std::shared_ptr> timeLoop, - const Communication& comm = Dune::MPIHelper::getCollectiveCommunication(), - const std::string& paramGroup = "") - : comm_(comm) - , timeLoop_(timeLoop) - , endIterMsgStream_(std::ostringstream::out) - , paramGroup_(paramGroup) - { - initParams_(paramGroup); - } - - //! the communicator for parallel runs - const Communication& comm() const - { return comm_; } - - /*! - * \brief Set the maximum acceptable difference of any primary variable - * between two iterations for declaring convergence. - * - * \param tolerance The maximum relative shift between two Newton - * iterations at which the scheme is considered finished - */ - void setMaxRelativeShift(Scalar tolerance) - { shiftTolerance_ = tolerance; } - - /*! - * \brief Set the maximum acceptable absolute residual for declaring convergence. - * - * \param tolerance The maximum absolute residual at which - * the scheme is considered finished - */ - void setMaxAbsoluteResidual(Scalar tolerance) - { residualTolerance_ = tolerance; } - - /*! - * \brief Set the maximum acceptable residual norm reduction. - * - * \param tolerance The maximum reduction of the residual norm - * at which the scheme is considered finished - */ - void setResidualReduction(Scalar tolerance) - { reductionTolerance_ = tolerance; } - - /*! - * \brief Set the number of iterations at which the Newton method - * should aim at. - * - * This is used to control the time-step size. The heuristic used - * is to scale the last time-step size by the deviation of the - * number of iterations used from the target steps. - * - * \param targetSteps Number of iterations which are considered "optimal" - */ - void setTargetSteps(int targetSteps) - { targetSteps_ = targetSteps; } - - /*! - * \brief Set the number of iterations after which the Newton - * method gives up. - * - * \param maxSteps Number of iterations after we give up - */ - void setMaxSteps(int maxSteps) - { maxSteps_ = maxSteps; } - - /*! - * \brief Returns true if another iteration should be done. - * - * \param uCurrentIter The solution of the current Newton iteration - * \param converged if the Newton method's convergence criterion was met in this step - */ - template - bool newtonProceed(const SolutionVector &uCurrentIter, bool converged) - { - if (numSteps_ < 2) - return true; // we always do at least two iterations - else if (converged) { - return false; // we are below the desired tolerance - } - else if (numSteps_ >= maxSteps_) { - // We have exceeded the allowed number of steps. If the - // maximum relative shift was reduced by a factor of at least 4, - // we proceed even if we are above the maximum number of steps. - if (enableShiftCriterion_) - return shift_*4.0 < lastShift_; - else - return reduction_*4.0 < lastReduction_; - } - - return true; - } - - /*! - * \brief Returns true if the error of the solution is below the - * tolerance. - */ - bool newtonConverged() const - { - if (enableShiftCriterion_ && !enableResidualCriterion_) - { - return shift_ <= shiftTolerance_; - } - else if (!enableShiftCriterion_ && enableResidualCriterion_) - { - if(enableAbsoluteResidualCriterion_) - return residualNorm_ <= residualTolerance_; - else - return reduction_ <= reductionTolerance_; - } - else if (satisfyResidualAndShiftCriterion_) - { - if(enableAbsoluteResidualCriterion_) - return shift_ <= shiftTolerance_ - && residualNorm_ <= residualTolerance_; - else - return shift_ <= shiftTolerance_ - && reduction_ <= reductionTolerance_; - } - else - { - return shift_ <= shiftTolerance_ - || reduction_ <= reductionTolerance_ - || residualNorm_ <= residualTolerance_; - } - - return false; - } - - /*! - * \brief Called before the Newton method is applied to an - * non-linear system of equations. - * - * \param u The initial solution - */ - template - void newtonBegin(const SolutionVector& u) - { - numSteps_ = 0; - } - - /*! - * \brief Indicates the beginning of a Newton iteration. - */ - template - void newtonBeginStep(const SolutionVector& u) - { - lastShift_ = shift_; - if (numSteps_ == 0) - { - lastReduction_ = 1.0; - } - else - { - lastReduction_ = reduction_; - } - } - - /*! - * \brief Returns the number of steps done since newtonBegin() was - * called. - */ - int newtonNumSteps() - { return numSteps_; } - - /*! - * \brief Update the maximum relative shift of the solution compared to - * the previous iteration. - * - * \param uLastIter The current iterative solution - * \param deltaU The difference between the current and the next solution - */ - template - void newtonUpdateShift(const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - shift_ = 0; - - for (int i = 0; i < int(uLastIter.size()); ++i) { - typename SolutionVector::block_type uNewI = uLastIter[i]; - uNewI -= deltaU[i]; - - Scalar shiftAtDof = relativeShiftAtDof_(uLastIter[i], uNewI); - using std::max; - shift_ = max(shift_, shiftAtDof); - } - - if (comm().size() > 1) - shift_ = comm().max(shift_); - } - - /*! - * \brief Assemble the linear system of equations \f$\mathbf{A}x - b = 0\f$. - * - * \param assembler The jacobian assembler - * \param uCurrentIter The current iteration's solution vector - */ - template - void assembleLinearSystem(JacobianAssembler& assembler, - const SolutionVector& uCurrentIter) - { - assembler.assembleJacobianAndResidual(uCurrentIter); - } - - /*! - * \brief Solve the linear system of equations \f$\mathbf{A}x - b = 0\f$. - * - * \throws Dumux::NumericalProblem if the linear solver didn't - * converge. - * - * \param ls The linear solver to be used - * \param A The matrix of the linear system of equations - * \param x The vector which solves the linear system - * \param b The right hand side of the linear system - */ - template - void solveLinearSystem(LinearSolver& ls, - JacobianMatrix& A, - SolutionVector& x, - SolutionVector& b) - { - try { - if (numSteps_ == 0) - { - Scalar norm2 = b.two_norm2(); - if (comm().size() > 1) - norm2 = comm().sum(norm2); - - using std::sqrt; - initialResidual_ = sqrt(norm2); - } - - //! Copy into a standard block vector. - //! This is necessary for all model _not_ using a FieldVector as - //! primary variables vector in combination with UMFPack or SuperLU as their interfaces are hard coded - //! 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; - using BlockType = Dune::FieldVector; - Dune::BlockVector xTmp; xTmp.resize(b.size()); - Dune::BlockVector bTmp(xTmp); - for (unsigned int i = 0; i < b.size(); ++i) - for (unsigned int j = 0; j < blockSize; ++j) - bTmp[i][j] = b[i][j]; - - int converged = ls.solve(A, xTmp, bTmp); - - for (unsigned int i = 0; i < x.size(); ++i) - for (unsigned int j = 0; j < blockSize; ++j) - x[i][j] = xTmp[i][j]; - - // make sure all processes converged - int convergedRemote = converged; - if (comm().size() > 1) - convergedRemote = comm().min(converged); - - if (!converged) { - DUNE_THROW(NumericalProblem, - "Linear solver did not converge"); - } - else if (!convergedRemote) { - DUNE_THROW(NumericalProblem, - "Linear solver did not converge on a remote process"); - } - } - catch (const Dune::Exception &e) { - // make sure all processes converged - int converged = 0; - if (comm().size() > 1) - converged = comm().min(converged); - - NumericalProblem p; - p.message(e.what()); - throw p; - } - } - - /*! - * \brief Update the current solution with a delta vector. - * - * The error estimates required for the newtonConverged() and - * newtonProceed() methods should be updated inside this method. - * - * Different update strategies, such as line search and chopped - * updates can be implemented. The default behavior is just to - * subtract deltaU from uLastIter, i.e. - * \f[ u^{k+1} = u^k - \Delta u^k \f] - * - * \param assembler The assembler (needed for global residual evaluation) - * \param uCurrentIter The solution vector after the current iteration - * \param uLastIter The solution vector after the last iteration - * \param deltaU The delta as calculated from solving the linear - * system of equations. This parameter also stores - * the updated solution. - */ - template - void newtonUpdate(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - if (enableShiftCriterion_) - newtonUpdateShift(uLastIter, deltaU); - - if (useLineSearch_) - { - lineSearchUpdate_(assembler, uCurrentIter, uLastIter, deltaU); - } - else { - for (unsigned int i = 0; i < uLastIter.size(); ++i) { - uCurrentIter[i] = uLastIter[i]; - uCurrentIter[i] -= deltaU[i]; - } - - if (enableResidualCriterion_) - { - residualNorm_ = assembler.residualNorm(uCurrentIter); - reduction_ = residualNorm_; - reduction_ /= initialResidual_; - } - else - { - // If we get here, the convergence criterion does not require - // additional residual evalutions. Thus, the grid variables have - // not yet been updated to the new uCurrentIter. - assembler.gridVariables().update(uCurrentIter); - } - } - } - - /*! - * \brief Indicates that one Newton iteration was finished. - * - * \param assembler The jacobian assembler - * \param uCurrentIter The solution after the current Newton iteration - * \param uLastIter The solution at the beginning of the current Newton iteration - */ - template - void newtonEndStep(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter) - { - ++numSteps_; - - if (verbose()) - { - std::cout << "\rNewton iteration " << numSteps_ << " done"; - if (enableShiftCriterion_) - std::cout << ", maximum relative shift = " << shift_; - if (enableResidualCriterion_ && enableAbsoluteResidualCriterion_) - std::cout << ", residual = " << residualNorm_; - else if (enableResidualCriterion_) - std::cout << ", residual reduction = " << reduction_; - std::cout << endIterMsg().str() << "\n"; - } - endIterMsgStream_.str(""); - - // When the Newton iterations are done: ask the model to check whether it makes sense - // TODO: how do we realize this? -> do this here in the newton controller - // model_().checkPlausibility(); - } - - /*! - * \brief Called if the Newton method ended - * (not known yet if we failed or succeeded) - */ - void newtonEnd() {} - - /*! - * \brief Called if the Newton method ended succcessfully - * This method is called _after_ newtonEnd() - */ - void newtonSucceed() {} - - /*! - * \brief Called if the Newton method broke down. - * This method is called _after_ newtonEnd() - */ - template - void newtonFail(Assembler& assembler, SolutionVector& u) - { - if (!assembler.isStationaryProblem()) - { - // set solution to previous solution - u = assembler.prevSol(); - - // reset the grid variables to the previous solution - assembler.gridVariables().resetTimeStep(u); - - if (verbose()) - { - std::cout << "Newton solver did not converge with dt = " - << timeLoop_->timeStepSize() << " seconds. Retrying with time step of " - << timeLoop_->timeStepSize()/2 << " seconds\n"; - } - - // try again with dt = dt/2 - timeLoop_->setTimeStepSize(timeLoop_->timeStepSize()/2); - } - else - DUNE_THROW(Dune::MathError, "Newton solver did not converge"); - } - - /*! - * \brief Suggest a new time-step size based on the old time-step - * size. - * - * The default behavior is to suggest the old time-step size - * scaled by the ratio between the target iterations and the - * iterations required to actually solve the last time-step. - */ - Scalar suggestTimeStepSize(Scalar oldTimeStep) const - { - // be aggressive reducing the time-step size but - // conservative when increasing it. the rationale is - // that we want to avoid failing in the next Newton - // iteration which would require another linearization - // of the problem. - if (numSteps_ > targetSteps_) { - Scalar percent = Scalar(numSteps_ - targetSteps_)/targetSteps_; - return oldTimeStep/(1.0 + percent); - } - - Scalar percent = Scalar(targetSteps_ - numSteps_)/targetSteps_; - return oldTimeStep*(1.0 + percent/1.2); - } - - std::ostringstream &endIterMsg() - { return endIterMsgStream_; } - - /*! - * \brief Specifies if the Newton method ought to be chatty. - */ - void setVerbose(bool val) - { verbose_ = val; } - - /*! - * \brief Returns true if the Newton method ought to be chatty. - */ - bool verbose() const - { return verbose_ && comm().rank() == 0; } - - /*! - * \brief Returns the parameter group - */ - const std::string& paramGroup() const - { return paramGroup_; } - -protected: - - //! initialize the parameters by reading from the parameter tree - void initParams_(const std::string& group = "") - { - useLineSearch_ = getParamFromGroup(group, "Newton.UseLineSearch"); - enableAbsoluteResidualCriterion_ = getParamFromGroup(group, "Newton.EnableAbsoluteResidualCriterion"); - enableShiftCriterion_ = getParamFromGroup(group, "Newton.EnableShiftCriterion"); - enableResidualCriterion_ = getParamFromGroup(group, "Newton.EnableResidualCriterion") || enableAbsoluteResidualCriterion_; - satisfyResidualAndShiftCriterion_ = getParamFromGroup(group, "Newton.SatisfyResidualAndShiftCriterion"); - if (!enableShiftCriterion_ && !enableResidualCriterion_) - { - DUNE_THROW(Dune::NotImplemented, - "at least one of NewtonEnableShiftCriterion or " - << "NewtonEnableResidualCriterion has to be set to true"); - } - - setMaxRelativeShift(getParamFromGroup(group, "Newton.MaxRelativeShift")); - setMaxAbsoluteResidual(getParamFromGroup(group, "Newton.MaxAbsoluteResidual")); - setResidualReduction(getParamFromGroup(group, "Newton.ResidualReduction")); - setTargetSteps(getParamFromGroup(group, "Newton.TargetSteps")); - setMaxSteps(getParamFromGroup(group, "Newton.MaxSteps")); - - verbose_ = true; - numSteps_ = 0; - } - - template - void lineSearchUpdate_(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - Scalar lambda = 1.0; - SolutionVector tmp(uLastIter); - - while (true) - { - uCurrentIter = deltaU; - uCurrentIter *= -lambda; - uCurrentIter += uLastIter; - - residualNorm_ = assembler.residualNorm(uCurrentIter); - reduction_ = residualNorm_; - reduction_ /= initialResidual_; - - if (reduction_ < lastReduction_ || lambda <= 0.125) { - endIterMsg() << ", residual reduction " << lastReduction_ << "->" << reduction_ << "@lambda=" << lambda; - return; - } - - // try with a smaller update - lambda /= 2.0; - } - } - - /*! - * \brief Returns the maximum relative shift between two vectors of - * primary variables. - * - * \param priVars1 The first vector of primary variables - * \param priVars2 The second vector of primary variables - */ - template - Scalar relativeShiftAtDof_(const PrimaryVariables &priVars1, - const PrimaryVariables &priVars2) - { - Scalar result = 0.0; - using std::abs; - using std::max; - // iterate over all primary variables - for (int j = 0; j < PrimaryVariables::dimension; ++j) { - Scalar eqErr = abs(priVars1[j] - priVars2[j]); - eqErr /= max(1.0,abs(priVars1[j] + priVars2[j])/2); - - result = max(result, eqErr); - } - return result; - } - - //! The communication object - Communication comm_; - - //! The time loop for stationary simulations - std::shared_ptr> timeLoop_; - - //! message stream to be displayed at the end of iterations - std::ostringstream endIterMsgStream_; - - //! switches on/off verbosity - bool verbose_; - - // shift criterion variables - Scalar shift_; - Scalar lastShift_; - Scalar shiftTolerance_; - - // residual criterion variables - Scalar reduction_; - Scalar residualNorm_; - Scalar lastReduction_; - Scalar initialResidual_; - Scalar reductionTolerance_; - Scalar residualTolerance_; - - //! optimal number of iterations we want to achieve - int targetSteps_; - //! maximum number of iterations we do before giving up - int maxSteps_; - //! actual number of steps done so far - int numSteps_; - - // further parameters - bool enablePartialReassemble_; - bool useLineSearch_; - bool enableAbsoluteResidualCriterion_; - bool enableShiftCriterion_; - bool enableResidualCriterion_; - bool satisfyResidualAndShiftCriterion_; - - //! the parameter group for getting parameters from the parameter tree - std::string paramGroup_; -}; - -} // end namespace Dumux - -#endif diff --git a/dumux/nonlinear/newtonmethod.hh b/dumux/nonlinear/newtonmethod.hh deleted file mode 100644 index 7264b3cedf..0000000000 --- a/dumux/nonlinear/newtonmethod.hh +++ /dev/null @@ -1,222 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/***************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup Nonlinear - * \brief The algorithmic part of the multi dimensional newton method. - * - * In order to use the method you need a Newtoncontroller - */ -#ifndef DUMUX_NEWTONMETHOD_HH -#define DUMUX_NEWTONMETHOD_HH - -#include -#include - -#include -#include - -#include -#include -#include -#include - -#warning "This file is deprecated. Use NewtonSolver instead." - -namespace Dumux { - -/*! - * \ingroup Newton - * \ingroup Nonlinear - * \brief The algorithmic part of the multi dimensional newton method. - * - * \tparam NewtonController the controller class is the driver implementation - * \tparam JacobianAssembler an assembler for the Jacobian and the residual - * \tparam LinearSolver the linear solver used to solve one iteration - */ -template -class DUNE_DEPRECATED_MSG("Use NewtonSolver instead.") NewtonMethod -{ - //! provide an interface as a form of type erasure - //! this is the minimal requirements a convergence write passed to a newton method has to fulfill - struct ConvergenceWriterInferface - { - template - void write(const SolutionVector &uLastIter, const SolutionVector &deltaU, const SolutionVector &residual) {} - }; - -public: - NewtonMethod(std::shared_ptr controller, - std::shared_ptr assembler, - std::shared_ptr linearSolver, - const std::string& modelParamGroup = "") - : controller_(controller) - , assembler_(assembler) - , linearSolver_(linearSolver) - { - // set the linear system (matrix & residual) in the assembler - assembler_->setLinearSystem(); - - // set a different default for the linear solver residual reduction - // within the Newton the linear solver doesn't need to solve too exact - using Scalar = typename LinearSolver::Scalar; - linearSolver_->setResidualReduction(getParamFromGroup(modelParamGroup, "LinearSolver.ResidualReduction", 1e-6)); - } - - /*! - * \brief Run the newton method to solve a non-linear system. - * The controller is responsible for all the strategic decisions. - */ - template - bool solve(SolutionVector& uCurrentIter, const std::unique_ptr& convWriter = nullptr) - { - try - { - // the given solution is the initial guess - SolutionVector uLastIter(uCurrentIter); - SolutionVector deltaU(uCurrentIter); - - Dune::Timer assembleTimer(false); - Dune::Timer solveTimer(false); - Dune::Timer updateTimer(false); - - // tell the controller that we begin solving - controller_->newtonBegin(uCurrentIter); - - // execute the method as long as the controller thinks - // that we should do another iteration - while (controller_->newtonProceed(uCurrentIter, controller_->newtonConverged())) - { - // notify the controller that we're about to start - // a new timestep - controller_->newtonBeginStep(uCurrentIter); - - // make the current solution to the old one - if (controller_->newtonNumSteps() > 0) - uLastIter = uCurrentIter; - - if (controller_->verbose()) { - std::cout << "Assemble: r(x^k) = dS/dt + div F - q; M = grad r"; - std::cout.flush(); - } - - /////////////// - // assemble - /////////////// - - // linearize the problem at the current solution - assembleTimer.start(); - controller_->assembleLinearSystem(*assembler_, uCurrentIter); - assembleTimer.stop(); - - /////////////// - // linear solve - /////////////// - - // Clear the current line using an ansi escape - // sequence. for an explanation see - // http://en.wikipedia.org/wiki/ANSI_escape_code - const char clearRemainingLine[] = { 0x1b, '[', 'K', 0 }; - - if (controller_->verbose()) { - std::cout << "\rSolve: M deltax^k = r"; - std::cout << clearRemainingLine; - std::cout.flush(); - } - - // solve the resulting linear equation system - solveTimer.start(); - - // set the delta vector to zero before solving the linear system! - deltaU = 0; - // ask the controller to solve the linearized system - controller_->solveLinearSystem(*linearSolver_, - assembler_->jacobian(), - deltaU, - assembler_->residual()); - solveTimer.stop(); - - /////////////// - // update - /////////////// - if (controller_->verbose()) { - std::cout << "\rUpdate: x^(k+1) = x^k - deltax^k"; - std::cout << clearRemainingLine; - std::cout.flush(); - } - - updateTimer.start(); - // update the current solution (i.e. uOld) with the delta - // (i.e. u). The result is stored in u - controller_->newtonUpdate(*assembler_, uCurrentIter, uLastIter, deltaU); - updateTimer.stop(); - - // tell the controller that we're done with this iteration - controller_->newtonEndStep(*assembler_, uCurrentIter, uLastIter); - - // if a convergence writer was specified compute residual and write output - if (convWriter) - { - assembler_->assembleResidual(uCurrentIter); - convWriter->write(uLastIter, deltaU, assembler_->residual()); - } - } - - // tell controller we are done - controller_->newtonEnd(); - - // reset state if newton failed - if (!controller_->newtonConverged()) - { - controller_->newtonFail(*assembler_, uCurrentIter); - return false; - } - - // tell controller we converged successfully - controller_->newtonSucceed(); - - if (controller_->verbose()) { - const auto elapsedTot = assembleTimer.elapsed() + solveTimer.elapsed() + updateTimer.elapsed(); - std::cout << "Assemble/solve/update time: " - << assembleTimer.elapsed() << "(" << 100*assembleTimer.elapsed()/elapsedTot << "%)/" - << solveTimer.elapsed() << "(" << 100*solveTimer.elapsed()/elapsedTot << "%)/" - << updateTimer.elapsed() << "(" << 100*updateTimer.elapsed()/elapsedTot << "%)" - << "\n"; - } - return true; - - } - catch (const NumericalProblem &e) - { - if (controller_->verbose()) - std::cout << "Newton: Caught exception: \"" << e.what() << "\"\n"; - controller_->newtonFail(*assembler_, uCurrentIter); - return false; - } - } - -private: - std::shared_ptr controller_; - std::shared_ptr assembler_; - std::shared_ptr linearSolver_; -}; - -} // end namespace Dumux - -#endif diff --git a/dumux/nonlinear/staggerednewtoncontroller.hh b/dumux/nonlinear/staggerednewtoncontroller.hh deleted file mode 100644 index bedcf6d4ad..0000000000 --- a/dumux/nonlinear/staggerednewtoncontroller.hh +++ /dev/null @@ -1,299 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/***************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup Nonlinear - * \ingroup StaggeredDiscretization - * \brief A newton controller for staggered finite volume schemes - */ -#ifndef DUMUX_STAGGERED_NEWTON_CONTROLLER_HH -#define DUMUX_STAGGERED_NEWTON_CONTROLLER_HH - -#include -#include -#include -#include - -#include -#include -#include -#include -#include - -#warning "This file is deprecated. Use NewtonSolver instead." - -namespace Dumux { - -/*! - * \ingroup Nonlinear - * \ingroup StaggeredDiscretization - * \brief A newton controller for staggered finite volume schemes - */ - -template > -class DUNE_DEPRECATED_MSG("Use NewtonSolver instead.") -StaggeredNewtonController : public NewtonController -{ - using ParentType = NewtonController; - -public: - using ParentType::ParentType; - - /*! - * \brief Solve the linear system of equations \f$\mathbf{A}x - b = 0\f$. - * - * Throws Dumux::NumericalProblem if the linear solver didn't - * converge. - * - * If the linear solver doesn't accept multitype matrices we copy the matrix - * into a 1x1 block BCRS matrix for solving. - * - * \param ls the linear solver - * \param A The matrix of the linear system of equations - * \param x The vector which solves the linear system - * \param b The right hand side of the linear system - */ - template - void solveLinearSystem(LinearSolver& ls, - JacobianMatrix& A, - SolutionVector& x, - SolutionVector& b) - { - // check matrix sizes - assert(checkMatrix_(A) && "Sub blocks of MultiType matrix have wrong sizes!"); - - try - { - if (this->numSteps_ == 0) - { - Scalar norm2 = b.two_norm2(); - if (this->comm().size() > 1) - norm2 = this->comm().sum(norm2); - - using std::sqrt; - this->initialResidual_ = sqrt(norm2); - } - - // solve by calling the appropriate implementation depending on whether the linear solver - // is capable of handling MultiType matrices or not - bool converged = solveLinearSystem_(ls, A, x, b, - std::integral_constant::value>()); - - // make sure all processes converged - int convergedRemote = converged; - if (this->comm().size() > 1) - convergedRemote = this->comm().min(converged); - - if (!converged) { - DUNE_THROW(NumericalProblem, - "Linear solver did not converge"); - } - else if (!convergedRemote) { - DUNE_THROW(NumericalProblem, - "Linear solver did not converge on a remote process"); - } - } - catch (const Dune::Exception &e) { - // make sure all processes converged - int converged = 0; - if (this->comm().size() > 1) - converged = this->comm().min(converged); - - NumericalProblem p; - p.message(e.what()); - throw p; - } - } - - - - /*! - * \brief Update the current solution with a delta vector. - * - * The error estimates required for the newtonConverged() and - * newtonProceed() methods should be updated inside this method. - * - * Different update strategies, such as line search and chopped - * updates can be implemented. The default behavior is just to - * subtract deltaU from uLastIter, i.e. - * \f[ u^{k+1} = u^k - \Delta u^k \f] - * - * \param assembler The assembler for Jacobian and residual - * \param uCurrentIter The solution vector after the current iteration - * \param uLastIter The solution vector after the last iteration - * \param deltaU The delta as calculated from solving the linear - * system of equations. This parameter also stores - * the updated solution. - */ - template - void newtonUpdate(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - if (this->enableShiftCriterion_) - this->newtonUpdateShift(uLastIter, deltaU); - - if (this->useLineSearch_) - { - this->lineSearchUpdate_(assembler, uCurrentIter, uLastIter, deltaU); - } - else { - using namespace Dune::Hybrid; - forEach(integralRange(Dune::Hybrid::size(uLastIter)), [&](const auto dofTypeIdx) - { - for (unsigned int i = 0; i < uLastIter[dofTypeIdx].size(); ++i) - { - uCurrentIter[dofTypeIdx][i] = uLastIter[dofTypeIdx][i]; - uCurrentIter[dofTypeIdx][i] -= deltaU[dofTypeIdx][i]; - } - }); - - if (this->enableResidualCriterion_) - { - this->residualNorm_ = assembler.residualNorm(uCurrentIter); - this->reduction_ = this->residualNorm_; - this->reduction_ /= this->initialResidual_; - } - } - - // update the variables class to the new solution - assembler.gridVariables().update(uCurrentIter); - } - - /*! - * \brief Update the maximum relative shift of the solution compared to - * the previous iteration. - * - * \param uLastIter The current iterative solution - * \param deltaU The difference between the current and the next solution - */ - template - void newtonUpdateShift(const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - this->shift_ = 0; - - using namespace Dune::Hybrid; - forEach(integralRange(Dune::Hybrid::size(uLastIter)), [&](const auto dofTypeIdx) - { - for (int i = 0; i < int(uLastIter[dofTypeIdx].size()); ++i) - { - auto uNewI = uLastIter[dofTypeIdx][i]; - uNewI -= deltaU[dofTypeIdx][i]; - - Scalar shiftAtDof = this->relativeShiftAtDof_(uLastIter[dofTypeIdx][i], uNewI); - this->shift_ = std::max(this->shift_, shiftAtDof); - } - }); - - if (this->comm().size() > 1) - this->shift_ = this->comm().max(this->shift_); - } - -private: - /*! - * \brief Solve the linear system of equations \f$\mathbf{A}x - b = 0\f$. - * - * Throws Dumux::NumericalProblem if the linear solver didn't - * converge. - * - * Specialization for linear solvers that can handle MultiType matrices. - * - */ - template - bool solveLinearSystem_(LinearSolver& ls, - JacobianMatrix& A, - SolutionVector& x, - SolutionVector& b, - std::true_type) - { - // TODO: automatically derive the precondBlockLevel - return ls.template solve(A, x, b); - } - - /*! - * \brief Solve the linear system of equations \f$\mathbf{A}x - b = 0\f$. - * - * Throws Dumux::NumericalProblem if the linear solver didn't - * converge. - * - * Specialization for linear solvers that cannot handle MultiType matrices. - * We copy the matrix into a 1x1 block BCRS matrix before solving. - * - */ - template - bool solveLinearSystem_(LinearSolver& ls, - JacobianMatrix& A, - SolutionVector& x, - SolutionVector& b, - std::false_type) - { - // create the bcrs matrix the IterativeSolver backend can handle - const auto M = MatrixConverter::multiTypeToBCRSMatrix(A); - - // get the new matrix sizes - const std::size_t numRows = M.N(); - assert(numRows == M.M()); - - // create the vector the IterativeSolver backend can handle - const auto bTmp = VectorConverter::multiTypeToBlockVector(b); - assert(bTmp.size() == numRows); - - // create a blockvector to which the linear solver writes the solution - using VectorBlock = typename Dune::FieldVector; - using BlockVector = typename Dune::BlockVector; - BlockVector y(numRows); - - // solve - const bool converged = ls.solve(M, y, bTmp); - - // copy back the result y into x - if(converged) - VectorConverter::retrieveValues(x, y); - - return converged; - } - - //! helper method to assure the MultiType matrix's sub blocks have the correct sizes - template - bool checkMatrix_(const JacobianMatrix& A) - { - bool matrixHasCorrectSize = true; - using namespace Dune::Hybrid; - using namespace Dune::Indices; - forEach(A, [&matrixHasCorrectSize](const auto& rowOfMultiTypeMatrix) - { - const auto numRowsLeftMostBlock = rowOfMultiTypeMatrix[_0].N(); - - forEach(rowOfMultiTypeMatrix, [&matrixHasCorrectSize, &numRowsLeftMostBlock](const auto& subBlock) - { - if (subBlock.N() != numRowsLeftMostBlock) - matrixHasCorrectSize = false; - }); - }); - return matrixHasCorrectSize; - } - -}; - -} // end namespace Dumux - -#endif diff --git a/dumux/porousmediumflow/2pnc/CMakeLists.txt b/dumux/porousmediumflow/2pnc/CMakeLists.txt index 612765aa72..5e3c930c89 100644 --- a/dumux/porousmediumflow/2pnc/CMakeLists.txt +++ b/dumux/porousmediumflow/2pnc/CMakeLists.txt @@ -3,7 +3,6 @@ install(FILES indices.hh model.hh -newtoncontroller.hh primaryvariableswitch.hh volumevariables.hh vtkoutputfields.hh diff --git a/dumux/porousmediumflow/compositional/privarswitchnewtoncontroller.hh b/dumux/porousmediumflow/compositional/privarswitchnewtoncontroller.hh deleted file mode 100644 index 553414dbbc..0000000000 --- a/dumux/porousmediumflow/compositional/privarswitchnewtoncontroller.hh +++ /dev/null @@ -1,220 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/**************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup PorousmediumCompositional - * \brief Reference implementation of a controller class for the Newton solver. - * - * Usually this controller should be sufficient. - */ -#ifndef DUMUX_PRIVARSWITCH_NEWTON_CONTROLLER_HH -#define DUMUX_PRIVARSWITCH_NEWTON_CONTROLLER_HH - -#include - -#include -#include -#include -#include -#include -#include - -#warning "This file is deprecated. Use PriVarSwitchNewtonSolver instead." - -namespace Dumux { - -/*! - * \ingroup PorousmediumCompositional - * \brief A newton controller that handles primary variable switches - * \todo make this independent of TypeTag by making PrimaryVariableSwitch a template argument - * and extracting everything model specific from there - * \todo Implement for volume variable caching enabled - */ -template -class DUNE_DEPRECATED_MSG("Use PriVarSwitchNewtonSolver instead.") -PriVarSwitchNewtonController : public NewtonController -{ - using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); - using ParentType = NewtonController; - using PrimaryVariableSwitch = typename GET_PROP_TYPE(TypeTag, PrimaryVariableSwitch); - - static constexpr bool isBox = GET_PROP_TYPE(TypeTag, FVGridGeometry)::discMethod == DiscretizationMethod::box; - -public: - using ParentType::ParentType; - - /*! - * \brief Returns true if the error of the solution is below the - * tolerance. - */ - bool newtonConverged() const - { - if (switchedInLastIteration_) - return false; - - return ParentType::newtonConverged(); - } - - /*! - * - * \brief Called before the Newton method is applied to an - * non-linear system of equations. - * - * \param u The initial solution - */ - template - void newtonBegin(const SolutionVector &u) - { - ParentType::newtonBegin(u); - priVarSwitch_ = std::make_unique(u.size()); - } - - /*! - * \brief Indicates that one Newton iteration was finished. - * - * \param assembler The jacobian assembler - * \param uCurrentIter The solution after the current Newton iteration - * \param uLastIter The solution at the beginning of the current Newton iteration - */ - template - void newtonEndStep(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter) - { - ParentType::newtonEndStep(assembler, uCurrentIter, uLastIter); - - // update the variable switch (returns true if the pri vars at at least one dof were switched) - // for disabled grid variable caching - const auto& fvGridGeometry = assembler.fvGridGeometry(); - const auto& problem = assembler.problem(); - auto& gridVariables = assembler.gridVariables(); - - // invoke the primary variable switch - switchedInLastIteration_ = priVarSwitch_->update(uCurrentIter, gridVariables, - problem, fvGridGeometry); - - if(switchedInLastIteration_) - { - for (const auto& element : elements(fvGridGeometry.gridView())) - { - // if the volume variables are cached globally, we need to update those where the primary variables have been switched - updateSwitchedVolVars_(std::integral_constant(), - element, assembler, uCurrentIter, uLastIter); - - // if the flux variables are cached globally, we need to update those where the primary variables have been switched - // (not needed for box discretization) - updateSwitchedFluxVarsCache_(std::integral_constant(), - element, assembler, uCurrentIter, uLastIter); - } - } - } - - /*! - * \brief Called if the Newton method ended - * (not known yet if we failed or succeeded) - */ - void newtonEnd() - { - ParentType::newtonEnd(); - - // in any way, we have to reset the switch flag - switchedInLastIteration_ = false; - // free some memory - priVarSwitch_.release(); - } - -private: - - /*! - * \brief Update the volume variables whose primary variables were - switched. Required when volume variables are cached globally. - */ - template - void updateSwitchedVolVars_(std::true_type, - const Element& element, - JacobianAssembler& assembler, - const SolutionVector &uCurrentIter, - const SolutionVector &uLastIter) - { - const auto& fvGridGeometry = assembler.fvGridGeometry(); - const auto& problem = assembler.problem(); - auto& gridVariables = assembler.gridVariables(); - - // make sure FVElementGeometry is bound to the element - auto fvGeometry = localView(fvGridGeometry); - fvGeometry.bindElement(element); - - // update the secondary variables if global caching is enabled - for (auto&& scv : scvs(fvGeometry)) - { - const auto dofIdxGlobal = scv.dofIndex(); - if (priVarSwitch_->wasSwitched(dofIdxGlobal)) - { - const auto elemSol = elementSolution(element, uCurrentIter, fvGridGeometry); - auto& volVars = gridVariables.curGridVolVars().volVars(scv); - volVars.update(elemSol, problem, element, scv); - } - } - } - - /*! - * \brief Update the fluxVars cache for dof whose primary variables were - switched. Required when flux variables are cached globally. - */ - template - void updateSwitchedFluxVarsCache_(std::true_type, - const Element& element, - JacobianAssembler& assembler, - const SolutionVector &uCurrentIter, - const SolutionVector &uLastIter) - { - const auto& fvGridGeometry = assembler.fvGridGeometry(); - auto& gridVariables = assembler.gridVariables(); - - // update the flux variables if global caching is enabled - const auto dofIdxGlobal = fvGridGeometry.dofMapper().index(element); - - if (priVarSwitch_->wasSwitched(dofIdxGlobal)) - { - // make sure FVElementGeometry and the volume variables are bound - auto fvGeometry = localView(fvGridGeometry); - fvGeometry.bind(element); - auto curElemVolVars = localView(gridVariables.curGridVolVars()); - curElemVolVars.bind(element, fvGeometry, uCurrentIter); - gridVariables.gridFluxVarsCache().updateElement(element, fvGeometry, curElemVolVars); - } - } - - //! brief Do nothing when volume variables are not cached globally. - template - void updateSwitchedVolVars_(std::false_type, Args&&... args) const {} - - //! brief Do nothing when flux variables are not cached globally. - template - void updateSwitchedFluxVarsCache_(std::false_type, Args&&... args) const {} - - //! the class handling the primary variable switch - std::unique_ptr priVarSwitch_; - //! if we switched primary variables in the last iteration - bool switchedInLastIteration_ = false; -}; - -} // end namespace Dumux - -#endif diff --git a/dumux/porousmediumflow/nonequilibrium/CMakeLists.txt b/dumux/porousmediumflow/nonequilibrium/CMakeLists.txt index 31fa8eb04f..09252e0484 100644 --- a/dumux/porousmediumflow/nonequilibrium/CMakeLists.txt +++ b/dumux/porousmediumflow/nonequilibrium/CMakeLists.txt @@ -4,8 +4,8 @@ install(FILES localresidual.hh gridvariables.hh indices.hh -newtoncontroller.hh volumevariables.hh vtkoutputfields.hh model.hh +newtonsolver.hh DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/porousmediumflow/nonequilibrium) diff --git a/dumux/porousmediumflow/nonequilibrium/newtoncontroller.hh b/dumux/porousmediumflow/nonequilibrium/newtoncontroller.hh deleted file mode 100644 index b1bfec283b..0000000000 --- a/dumux/porousmediumflow/nonequilibrium/newtoncontroller.hh +++ /dev/null @@ -1,66 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/***************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup PorousmediumNonEquilibriumModel - * \brief A MpNc specific controller for the newton solver. - * This controller calls the velocity averaging in the model after each iteration. - */ -#ifndef DUMUX_NONEQUILIBRIUM_NEWTON_CONTROLLER_HH -#define DUMUX_NONEQUILIBRIUM_NEWTON_CONTROLLER_HH - -#include - -#include - -namespace Dumux { -/*! - * \ingroup PorousmediumNonEquilibriumModel - * \brief A nonequilibrium specific controller for the newton solver. - * This controller calls the velocity averaging in the problem after each iteration. - */ -template > -class NonEquilibriumNewtonController : public NewtonController -{ - using ParentType = NewtonController; - -public: - using ParentType::ParentType; - - template - void newtonUpdate(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - ParentType::newtonUpdate(assembler, - uCurrentIter, - uLastIter, - deltaU); - - auto& gridVariables = assembler.gridVariables(); - // Averages the face velocities of a vertex. Implemented in the model. - // The velocities are stored in the model. - gridVariables.calcVelocityAverage(uCurrentIter); - } -}; - -} // end namespace Dumux -#endif // DUMUX_VELO_PROB_AVERAGE_NEWTON_CONTROLLER_HH diff --git a/dumux/porousmediumflow/richards/CMakeLists.txt b/dumux/porousmediumflow/richards/CMakeLists.txt index 2929d226a4..199fc1455f 100644 --- a/dumux/porousmediumflow/richards/CMakeLists.txt +++ b/dumux/porousmediumflow/richards/CMakeLists.txt @@ -4,7 +4,8 @@ install(FILES indices.hh localresidual.hh model.hh -newtoncontroller.hh +newtonsolver.hh +privarswitchnewtonsolver.hh primaryvariableswitch.hh volumevariables.hh vtkoutputfields.hh diff --git a/dumux/porousmediumflow/richards/newtoncontroller.hh b/dumux/porousmediumflow/richards/newtoncontroller.hh deleted file mode 100644 index e6a7cccd48..0000000000 --- a/dumux/porousmediumflow/richards/newtoncontroller.hh +++ /dev/null @@ -1,122 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/***************************************************************************** - * See the file COPYING for full copying permissions. * - * * - * This program is free software: you can redistribute it and/or modify * - * it under the terms of the GNU General Public License as published by * - * the Free Software Foundation, either version 2 of the License, or * - * (at your option) any later version. * - * * - * This program is distributed in the hope that it will be useful, * - * but WITHOUT ANY WARRANTY; without even the implied warranty of * - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * - * GNU General Public License for more details. * - * * - * You should have received a copy of the GNU General Public License * - * along with this program. If not, see . * - *****************************************************************************/ -/*! - * \file - * \ingroup RichardsModel - * \brief A Richards model specific controller for the newton solver. - */ -#ifndef DUMUX_RICHARDS_NEWTON_CONTROLLER_HH -#define DUMUX_RICHARDS_NEWTON_CONTROLLER_HH - -#include -#include -#include -#include - -#warning "This file is deprecated. Use RichardsNewtonSolver instead." - -namespace Dumux { -/*! - * \ingroup RichardsModel - * \brief A Richards model specific controller for the newton solver. - * - * This controller 'knows' what a 'physically meaningful' solution is - * and can thus do update smarter than the plain Newton controller. - * - * \todo make this typetag independent by extracting anything model specific from assembler - * or from possible ModelTraits. - */ -template -class DUNE_DEPRECATED_MSG("Use RichardsNewtonSolver instead.") -RichardsNewtonController : public NewtonController -{ - using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); - using ParentType = NewtonController; - - using MaterialLaw = typename GET_PROP_TYPE(TypeTag, MaterialLaw); - using Indices = typename GET_PROP_TYPE(TypeTag, Indices); - enum { pressureIdx = Indices::pressureIdx }; - -public: - using ParentType::ParentType; - - /*! - * \brief Update the current solution of the newton method - * - * This is basically the step - * \f[ u^{k+1} = u^k - \Delta u^k \f] - * - * \param assembler The Jacobian assembler - * \param uCurrentIter The solution after the current Newton iteration \f$ u^{k+1} \f$ - * \param uLastIter The solution after the last Newton iteration \f$ u^k \f$ - * \param deltaU The vector of differences between the last - * iterative solution and the next one \f$ \Delta u^k \f$ - */ - template - void newtonUpdate(JacobianAssembler& assembler, - SolutionVector &uCurrentIter, - const SolutionVector &uLastIter, - const SolutionVector &deltaU) - { - ParentType::newtonUpdate(assembler, uCurrentIter, uLastIter, deltaU); - if (!this->useLineSearch_ && getParamFromGroup(this->paramGroup(), "Newton.EnableChop")) - { - // do not clamp anything after 5 iterations - if (this->numSteps_ > 4) - return; - - // clamp saturation change to at most 20% per iteration - const auto& fvGridGeometry = assembler.fvGridGeometry(); - for (const auto& element : elements(fvGridGeometry.gridView())) - { - auto fvGeometry = localView(fvGridGeometry); - fvGeometry.bindElement(element); - - for (auto&& scv : scvs(fvGeometry)) - { - auto dofIdxGlobal = scv.dofIndex(); - - // calculate the old wetting phase saturation - const auto& spatialParams = assembler.problem().spatialParams(); - const auto elemSol = elementSolution(element, uCurrentIter, fvGridGeometry); - const auto& materialLawParams = spatialParams.materialLawParams(element, scv, elemSol); - const Scalar pcMin = MaterialLaw::pc(materialLawParams, 1.0); - const Scalar pw = uLastIter[dofIdxGlobal][pressureIdx]; - using std::max; - const Scalar pn = max(assembler.problem().nonWettingReferencePressure(), pw + pcMin); - const Scalar pcOld = pn - pw; - const Scalar SwOld = max(0.0, MaterialLaw::sw(materialLawParams, pcOld)); - - // convert into minimum and maximum wetting phase - // pressures - const Scalar pwMin = pn - MaterialLaw::pc(materialLawParams, SwOld - 0.2); - const Scalar pwMax = pn - MaterialLaw::pc(materialLawParams, SwOld + 0.2); - - // clamp the result - using std::min; using std::max; - uCurrentIter[dofIdxGlobal][pressureIdx] = max(pwMin, min(uCurrentIter[dofIdxGlobal][pressureIdx], pwMax)); - } - } - } - } -}; - -} // end namespace Dumux - -#endif diff --git a/dumux/porousmediumflow/richardsnc/model.hh b/dumux/porousmediumflow/richardsnc/model.hh index 659a42feb0..91720a3e1c 100644 --- a/dumux/porousmediumflow/richardsnc/model.hh +++ b/dumux/porousmediumflow/richardsnc/model.hh @@ -70,7 +70,6 @@ #include #include -#include #include #include -- GitLab From e56fc2e45a04a1e82e4a46b8d63bf39ebd146eff Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Mon, 9 Apr 2018 07:27:27 +0200 Subject: [PATCH 19/19] [NewtonSolver] Fix docu --- dumux/nonlinear/newtonsolver.hh | 27 ++++++++++--------- dumux/nonlinear/privarswitchnewtonsolver.hh | 8 ++---- .../nonequilibrium/newtonsolver.hh | 2 +- .../porousmediumflow/richards/newtonsolver.hh | 2 +- 4 files changed, 18 insertions(+), 21 deletions(-) diff --git a/dumux/nonlinear/newtonsolver.hh b/dumux/nonlinear/newtonsolver.hh index b50dbd7c63..82f9eca3de 100644 --- a/dumux/nonlinear/newtonsolver.hh +++ b/dumux/nonlinear/newtonsolver.hh @@ -19,9 +19,9 @@ /*! * \file * \ingroup Nonlinear - * \brief Reference implementation of a controller class for the Newton solver. + * \brief Reference implementation of the Newton solver. * - * Usually this controller should be sufficient. + * Usually this solver should be sufficient. */ #ifndef DUMUX_NEWTON_SOLVER_HH #define DUMUX_NEWTON_SOLVER_HH @@ -61,11 +61,12 @@ struct supportsPartialReassembly /*! * \ingroup Nonlinear - * \brief An implementation of a Newton controller - * \tparam Scalar the scalar type + * \brief An implementation of a Newton solver + * \tparam Assembler the assembler + * \tparam LinearSolver the linear solver * \tparam Comm the communication object used to communicate with all processes * \note If you want to specialize only some methods but are happy with the - * defaults of the reference controller, derive your controller from + * defaults of the reference solver, derive your solver from * this class and simply overload the required methods. */ template & convWriter = nullptr) { @@ -471,7 +472,7 @@ public: endIterMsgStream_.str(""); // When the Newton iterations are done: ask the model to check whether it makes sense - // TODO: how do we realize this? -> do this here in the newton controller + // TODO: how do we realize this? -> do this here in the newton solver // model_().checkPlausibility(); } @@ -620,7 +621,7 @@ private: /*! * \brief Run the newton method to solve a non-linear system. - * The controller is responsible for all the strategic decisions. + * The solver is responsible for all the strategic decisions. */ bool solve_(SolutionVector& uCurrentIter, const std::unique_ptr& convWriter = nullptr) { @@ -643,11 +644,11 @@ private: { newtonBegin(uCurrentIter); - // execute the method as long as the controller thinks + // execute the method as long as the solver thinks // that we should do another iteration while (newtonProceed(uCurrentIter, newtonConverged())) { - // notify the controller that we're about to start + // notify the solver that we're about to start // a new timestep newtonBeginStep(uCurrentIter); @@ -708,7 +709,7 @@ private: newtonUpdate(uCurrentIter, uLastIter, deltaU); updateTimer.stop(); - // tell the controller that we're done with this iteration + // tell the solver that we're done with this iteration newtonEndStep(uCurrentIter, uLastIter); // if a convergence writer was specified compute residual and write output @@ -719,7 +720,7 @@ private: } } - // tell controller we are done + // tell solver we are done newtonEnd(); // reset state if newton failed @@ -729,7 +730,7 @@ private: return false; } - // tell controller we converged successfully + // tell solver we converged successfully newtonSucceed(); if (verbose_) { diff --git a/dumux/nonlinear/privarswitchnewtonsolver.hh b/dumux/nonlinear/privarswitchnewtonsolver.hh index 39359fdc1e..f28bb46a38 100644 --- a/dumux/nonlinear/privarswitchnewtonsolver.hh +++ b/dumux/nonlinear/privarswitchnewtonsolver.hh @@ -19,9 +19,7 @@ /*! * \file * \ingroup PorousmediumCompositional - * \brief Reference implementation of a controller class for the Newton solver. - * - * Usually this controller should be sufficient. + * \copydoc Dumux::PriVarSwitchNewtonSolver */ #ifndef DUMUX_PRIVARSWITCH_NEWTON_SOLVER_HH #define DUMUX_PRIVARSWITCH_NEWTON_SOLVER_HH @@ -37,9 +35,7 @@ namespace Dumux { /*! * \ingroup PorousmediumCompositional - * \brief A newton controller that handles primary variable switches - * \todo make this independent of TypeTag by making PrimaryVariableSwitch a template argument - * and extracting everything model specific from there + * \brief A newton solver that handles primary variable switches */ template class PriVarSwitchNewtonSolver : public NewtonSolver diff --git a/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh b/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh index f12622d634..838d2f5e1c 100644 --- a/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh +++ b/dumux/porousmediumflow/nonequilibrium/newtonsolver.hh @@ -20,7 +20,7 @@ * \file * \ingroup PorousmediumNonEquilibriumModel * \brief A MpNc specific newton solver. - * This controller calls the velocity averaging in the model after each iteration. + * This solver calls the velocity averaging in the model after each iteration. */ #ifndef DUMUX_NONEQUILIBRIUM_NEWTON_SOLVER_HH #define DUMUX_NONEQUILIBRIUM_NEWTON_SOLVER_HH diff --git a/dumux/porousmediumflow/richards/newtonsolver.hh b/dumux/porousmediumflow/richards/newtonsolver.hh index 1cd9e43786..0c89333066 100644 --- a/dumux/porousmediumflow/richards/newtonsolver.hh +++ b/dumux/porousmediumflow/richards/newtonsolver.hh @@ -33,7 +33,7 @@ namespace Dumux { * \ingroup RichardsModel * \brief A Richards model specific newton solver. * - * This controller 'knows' what a 'physically meaningful' solution is + * This solver 'knows' what a 'physically meaningful' solution is * and can thus do update smarter than the plain Newton solver. * * \todo make this typetag independent by extracting anything model specific from assembler -- GitLab