Commit 2ffab8e9 authored by Timo Koch's avatar Timo Koch Committed by Dennis Gläser
Browse files

[newton] Extract privarswitch implementation to separate class

parent 0101e5a9
......@@ -2,5 +2,6 @@ install(FILES
findscalarroot.hh
newtonconvergencewriter.hh
newtonsolver.hh
primaryvariableswitchadapter.hh
staggerednewtonconvergencewriter.hh
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/nonlinear)
......@@ -55,10 +55,32 @@
#include "newtonconvergencewriter.hh"
#include "newtonvariablesbackend.hh"
#include "primaryvariableswitchadapter.hh"
namespace Dumux {
namespace Detail {
// Helper boolean that states if the assembler exports grid variables
template<class Assembler> using AssemblerGridVariablesType = typename Assembler::GridVariables;
template<class Assembler>
inline constexpr bool assemblerExportsGridVariables
= Dune::Std::is_detected_v<AssemblerGridVariablesType, Assembler>;
// helper struct to define the variables on which the privarswitch should operate
template<class Assembler, bool exportsGridVars = assemblerExportsGridVariables<Assembler>>
struct PriVarSwitchVariablesType { using Type = typename Assembler::GridVariables; };
// if assembler does not export them, use an empty class. These situations either mean
// that there is no privarswitch, or, it is handled by a derived implementation.
template<class Assembler>
struct PriVarSwitchVariablesType<Assembler, false>
{ using Type = struct EmptyGridVariables {}; };
// Helper alias to deduce the variables types used in the privarswitch adapter
template<class Assembler>
using PriVarSwitchVariables
= typename PriVarSwitchVariablesType<Assembler, assemblerExportsGridVariables<Assembler>>::Type;
//! helper struct detecting if an assembler supports partial reassembly
struct supportsPartialReassembly
{
......@@ -69,13 +91,6 @@ struct supportsPartialReassembly
{}
};
//! helper aliases to extract a primary variable switch from the VolumeVariables (if defined, yields int otherwise)
template<class Assembler>
using DetectPVSwitch = typename Assembler::GridVariables::VolumeVariables::PrimaryVariableSwitch;
template<class Assembler>
using GetPVSwitch = Dune::Std::detected_or<int, DetectPVSwitch, Assembler>;
// helper struct and function detecting if the linear solver features a norm() function
template <class LinearSolver, class Residual>
using NormDetector = decltype(std::declval<LinearSolver>().norm(std::declval<Residual>()));
......@@ -188,10 +203,14 @@ private:
using ConvergenceWriter = ConvergenceWriterInterface<SolutionVector>;
using TimeLoop = TimeLoopBase<Scalar>;
using PrimaryVariableSwitch = typename Detail::GetPVSwitch<Assembler>::type;
using HasPriVarsSwitch = typename Detail::GetPVSwitch<Assembler>::value_t; // std::true_type or std::false_type
static constexpr bool hasPriVarsSwitch() { return HasPriVarsSwitch::value; };
// enable models with primary variable switch
// TODO: Always use ParentType::Variables once we require assemblers to export variables
static constexpr bool assemblerExportsVariables = Impl::exportsVariables<Assembler>;
using PriVarSwitchVariables
= std::conditional_t<assemblerExportsVariables,
typename ParentType::Variables,
Detail::PriVarSwitchVariables<Assembler>>;
using PrimaryVariableSwitchAdapter = Dumux::PrimaryVariableSwitchAdapter<PriVarSwitchVariables>;
public:
using typename ParentType::Variables;
......@@ -208,6 +227,7 @@ public:
, endIterMsgStream_(std::ostringstream::out)
, comm_(comm)
, paramGroup_(paramGroup)
, priVarSwitchAdapter_(std::make_unique<PrimaryVariableSwitchAdapter>(paramGroup))
{
initParams_(paramGroup);
......@@ -221,12 +241,6 @@ public:
// initialize the partial reassembler
if (enablePartialReassembly_)
partialReassembler_ = std::make_unique<Reassembler>(this->assembler());
if (hasPriVarsSwitch())
{
const int priVarSwitchVerbosity = getParamFromGroup<int>(paramGroup, "PrimaryVariableSwitch.Verbosity", 1);
priVarSwitch_ = std::make_unique<PrimaryVariableSwitch>(priVarSwitchVerbosity);
}
}
//! the communicator for parallel runs
......@@ -370,13 +384,23 @@ public:
virtual void newtonBegin(Variables& initVars)
{
numSteps_ = 0;
initPriVarSwitch_(initVars, HasPriVarsSwitch{});
if constexpr (hasPriVarsSwitch<PriVarSwitchVariables>)
{
if constexpr (assemblerExportsVariables)
priVarSwitchAdapter_->initialize(Backend::getDofVector(initVars), initVars);
else // this assumes assembly with solution (i.e. Variables=SolutionVector)
priVarSwitchAdapter_->initialize(initVars, this->assembler().gridVariables());
}
const auto& initSol = Backend::getDofVector(initVars);
// write the initial residual if a convergence writer was set
if (convergenceWriter_)
{
this->assembler().assembleResidual(initVars);
// dummy vector, there is no delta before solving the linear system
auto delta = Backend::makeZeroDofVector(Backend::size(initSol));
convergenceWriter_->write(initSol, delta, this->assembler().residual());
}
......@@ -591,7 +615,13 @@ public:
virtual void newtonEndStep(Variables &vars,
const SolutionVector &uLastIter)
{
invokePriVarSwitch_(vars, HasPriVarsSwitch{});
if constexpr (hasPriVarsSwitch<PriVarSwitchVariables>)
{
if constexpr (assemblerExportsVariables)
priVarSwitchAdapter_->invoke(Backend::getDofVector(vars), vars);
else // this assumes assembly with solution (i.e. Variables=SolutionVector)
priVarSwitchAdapter_->invoke(vars, this->assembler().gridVariables());
}
++numSteps_;
......@@ -633,7 +663,7 @@ public:
{
// in case the model has a priVar switch and some some primary variables
// actually switched their state in the last iteration, enforce another iteration
if (priVarsSwitchedInLastIteration_)
if (priVarSwitchAdapter_->switched())
return false;
if (enableShiftCriterion_ && !enableResidualCriterion_)
......@@ -854,87 +884,6 @@ protected:
bool enableResidualCriterion() const
{ return enableResidualCriterion_; }
/*!
* \brief Initialize the privar switch, noop if there is no priVarSwitch
*/
void initPriVarSwitch_(Variables&, std::false_type) {}
/*!
* \brief Initialize the privar switch
*/
void initPriVarSwitch_(Variables& vars, std::true_type)
{
priVarSwitch_->reset(Backend::size(Backend::getDofVector(vars)));
priVarsSwitchedInLastIteration_ = false;
const auto& problem = this->assembler().problem();
const auto& gridGeometry = this->assembler().gridGeometry();
// if the assembler does not export the variables type
// we assume old-style assemblers which operate on solution
// vectors and return the grid variables
if constexpr (!assemblerExportsVariables)
priVarSwitch_->updateDirichletConstraints(problem, gridGeometry,
this->assembler().gridVariables(),
vars);
else
priVarSwitch_->updateDirichletConstraints(problem, gridGeometry, vars, Backend::getDofVector(vars));
}
/*!
* \brief Switch primary variables if necessary, noop if there is no priVarSwitch
*/
void invokePriVarSwitch_(Variables&, std::false_type) {}
/*!
* \brief Switch primary variables if necessary
*/
void invokePriVarSwitch_(Variables& vars, std::true_type)
{
// update the variable switch (returns true if the pri vars at at least one dof were switched)
// for disabled grid variable caching
const auto& gridGeometry = this->assembler().gridGeometry();
const auto& problem = this->assembler().problem();
// we assume an old-style assembler if it doesn't export the Variable type
if constexpr (!assemblerExportsVariables)
{
auto& gridVariables = this->assembler().gridVariables();
// invoke the primary variable switch
priVarsSwitchedInLastIteration_ = priVarSwitch_->update(vars, gridVariables, problem, gridGeometry);
if (priVarsSwitchedInLastIteration_)
{
for (const auto& element : elements(gridGeometry.gridView()))
{
// if the volume variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedVolVars(problem, element, gridGeometry, gridVariables, Backend::getDofVector(vars));
// if the flux variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedFluxVarsCache(problem, element, gridGeometry, gridVariables, Backend::getDofVector(vars));
}
}
}
else
{
// invoke the primary variable switch
priVarsSwitchedInLastIteration_ = priVarSwitch_->update(Backend::getDofVector(vars), vars, problem, gridGeometry);
if (priVarsSwitchedInLastIteration_)
{
for (const auto& element : elements(gridGeometry.gridView()))
{
// if the volume variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedVolVars(problem, element, gridGeometry, vars, Backend::getDofVector(vars));
// if the flux variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedFluxVarsCache(problem, element, gridGeometry, vars, Backend::getDofVector(vars));
}
}
}
}
//! optimal number of iterations we want to achieve
int targetSteps_;
//! minimum number of iterations we do
......@@ -1395,9 +1344,7 @@ private:
std::size_t numLinearSolverBreakdowns_ = 0; //! total number of linear solves that failed
//! the class handling the primary variable switch
std::unique_ptr<PrimaryVariableSwitch> priVarSwitch_;
//! if we switched primary variables in the last iteration
bool priVarsSwitchedInLastIteration_ = false;
std::unique_ptr<PrimaryVariableSwitchAdapter> priVarSwitchAdapter_;
//! convergence writer
std::shared_ptr<ConvergenceWriter> convergenceWriter_ = nullptr;
......
// -*- 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 3 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 <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup Nonlinear
* \brief An adapter for the Newton to manage models with primary variable switch
*/
#ifndef DUMUX_NONLINEAR_PRIMARY_VARIABLE_SWITCH_ADAPTER_HH
#define DUMUX_NONLINEAR_PRIMARY_VARIABLE_SWITCH_ADAPTER_HH
#include <memory>
#include <dune/common/std/type_traits.hh>
namespace Dumux {
namespace Detail {
//! helper aliases to extract a primary variable switch from the VolumeVariables (if defined, yields int otherwise)
template<class Variables>
using DetectPVSwitch = typename Variables::VolumeVariables::PrimaryVariableSwitch;
template<class Variables>
using GetPVSwitch = Dune::Std::detected_or<int, DetectPVSwitch, Variables>;
template<class Variables>
using PrimaryVariableSwitch = typename GetPVSwitch<Variables>::type;
} // end namespace Detail
/*!
* \ingroup Nonlinear
* \brief Helper boolean to check if the given variables involve primary variable switching.
*/
template<class Variables>
inline constexpr bool hasPriVarsSwitch = typename Detail::GetPVSwitch<Variables>::value_t();
/*!
* \ingroup Nonlinear
* \brief An adapter for the Newton to manage models with primary variable switch
*/
template <class Variables, bool isValid = hasPriVarsSwitch<Variables>>
class PrimaryVariableSwitchAdapter
{
using PrimaryVariableSwitch = typename Detail::PrimaryVariableSwitch<Variables>;
public:
PrimaryVariableSwitchAdapter(const std::string& paramGroup = "")
{
const int priVarSwitchVerbosity = getParamFromGroup<int>(paramGroup, "PrimaryVariableSwitch.Verbosity", 1);
priVarSwitch_ = std::make_unique<PrimaryVariableSwitch>(priVarSwitchVerbosity);
}
/*!
* \brief Initialize the privar switch
*/
template<class SolutionVector>
void initialize(SolutionVector& sol, Variables& vars)
{
priVarSwitch_->reset(sol.size());
priVarsSwitchedInLastIteration_ = false;
const auto& problem = vars.curGridVolVars().problem();
const auto& gridGeometry = problem.gridGeometry();
priVarSwitch_->updateDirichletConstraints(problem, gridGeometry, vars, sol);
}
/*!
* \brief Switch primary variables if necessary
*/
template<class SolutionVector>
void invoke(SolutionVector& uCurrentIter, Variables& vars)
{
// update the variable switch (returns true if the pri vars at at least one dof were switched)
// for disabled grid variable caching
const auto& problem = vars.curGridVolVars().problem();
const auto& gridGeometry = problem.gridGeometry();
// invoke the primary variable switch
priVarsSwitchedInLastIteration_ = priVarSwitch_->update(uCurrentIter, vars, problem, gridGeometry);
if (priVarsSwitchedInLastIteration_)
{
for (const auto& element : elements(gridGeometry.gridView()))
{
// if the volume variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedVolVars(problem, element, gridGeometry, vars, uCurrentIter);
// if the flux variables are cached globally, we need to update those where the primary variables have been switched
priVarSwitch_->updateSwitchedFluxVarsCache(problem, element, gridGeometry, vars, uCurrentIter);
}
}
}
/*!
* \brief Whether the primary variables have been switched in the last call to invoke
*/
bool switched() const
{ return priVarsSwitchedInLastIteration_; }
private:
//! the class handling the primary variable switch
std::unique_ptr<PrimaryVariableSwitch> priVarSwitch_;
//! if we switched primary variables in the last iteration
bool priVarsSwitchedInLastIteration_ = false;
};
/*!
* \ingroup Nonlinear
* \brief An empty adapter for the Newton for models without primary variable switch
*/
template <class Variables>
class PrimaryVariableSwitchAdapter<Variables, false>
{
public:
PrimaryVariableSwitchAdapter(const std::string& paramGroup = "") {}
template<class SolutionVector>
void initialize(SolutionVector&, Variables&) {}
template<class SolutionVector>
void invoke(SolutionVector&, Variables&) {}
bool switched() const { return false; }
};
} // end namespace Dumux
#endif
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment