diff --git a/dumux/multidomain/privarswitchnewtonsolver.hh b/dumux/multidomain/privarswitchnewtonsolver.hh new file mode 100644 index 0000000000000000000000000000000000000000..6397cb77aefc70d82120aca945c4c13e0307ae59 --- /dev/null +++ b/dumux/multidomain/privarswitchnewtonsolver.hh @@ -0,0 +1,160 @@ +// -*- 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 <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \ingroup MultiDomain + * \copydoc Dumux::MultiDomainPriVarSwitchNewtonSolver + */ +#ifndef DUMUX_MULTIDOMAIN_PRIVARSWITCH_NEWTON_SOLVER_HH +#define DUMUX_MULTIDOMAIN_PRIVARSWITCH_NEWTON_SOLVER_HH + +#include <memory> + +#include <dumux/common/typetraits/utility.hh> +#include <dumux/discretization/methods.hh> +#include <dumux/discretization/elementsolution.hh> +#include <dumux/multidomain/newtonsolver.hh> + +namespace Dumux { + +/*! + * \ingroup MultiDomain + * \brief A newton solver that handles primary variable switches for multi domain problems + */ +template <class Assembler, class LinearSolver, class CouplingManager, class PrivarSwitchTypeTuple, + class Reassembler = DefaultPartialReassembler, + class Comm = Dune::CollectiveCommunication<Dune::MPIHelper::MPICommunicator> > +class MultiDomainPriVarSwitchNewtonSolver : public MultiDomainNewtonSolver<Assembler, LinearSolver, CouplingManager, Reassembler, Comm> +{ + using Scalar = typename Assembler::Scalar; + using ParentType = MultiDomainNewtonSolver<Assembler, LinearSolver, CouplingManager, Reassembler, Comm>; + using SolutionVector = typename Assembler::ResidualType; + + using DomainIds = std::make_index_sequence<Assembler::Traits::numSubDomains>; + + template<std::size_t id> + using PrivarSwitchPtr = std::unique_ptr<std::tuple_element_t<id, PrivarSwitchTypeTuple>>; + + using PriVarSwitchPtrTuple = typename makeFromIndexedType<std::tuple, PrivarSwitchPtr, DomainIds>::type; + +public: + using ParentType::ParentType; + + /*! + * \brief Returns true if the error of the solution is below the + * tolerance. + */ + bool newtonConverged() const override + { + if (std::any_of(switchedInLastIteration_.begin(), switchedInLastIteration_.end(), [](bool i){ return i;})) + 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 + */ + void newtonBegin(const SolutionVector &u) override + { + ParentType::newtonBegin(u); + switchedInLastIteration_.fill(false); + + using namespace Dune::Hybrid; + forEach(DomainIds{}, [&](auto&& id) + { + std::get<id>(priVarSwitches_) = std::make_unique<std::tuple_element_t<id, PrivarSwitchTypeTuple>>(u[id].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 + */ + void newtonEndStep(SolutionVector &uCurrentIter, + const SolutionVector &uLastIter) override + { + ParentType::newtonEndStep(uCurrentIter, uLastIter); + + // update the variable switch (returns true if the pri vars at at least one dof were switched) + // for disabled grid variable caching + auto& assembler = this->assembler(); + + using namespace Dune::Hybrid; + forEach(DomainIds{}, [&](auto&& id) + { + const auto& fvGridGeometry = assembler.fvGridGeometry(id); + const auto& problem = assembler.problem(id); + auto& gridVariables = assembler.gridVariables(id); + auto& priVarSwitch = std::get<id>(priVarSwitches_); + + // invoke the primary variable switch + switchedInLastIteration_[id] = priVarSwitch->update(uCurrentIter[id], gridVariables, + problem, fvGridGeometry); + + if(switchedInLastIteration_[id]) + { + 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 + priVarSwitch->updateSwitchedVolVars(problem, element, fvGridGeometry, gridVariables, uCurrentIter); + + // if the flux variables are cached globally, we need to update those where the primary variables have been switched + priVarSwitch->updateSwitchedFluxVarsCache(problem, element, fvGridGeometry, gridVariables, uCurrentIter); + } + } + }); + } + + /*! + * \brief Called if the Newton method ended + * (not known yet if we failed or succeeded) + */ + void newtonEnd() override + { + ParentType::newtonEnd(); + + // in any way, we have to reset the switch flag + switchedInLastIteration_.fill(false); + + // free some memory + using namespace Dune::Hybrid; + forEach(DomainIds{}, [&](auto&& id) + { + std::get<id>(priVarSwitches_).release(); + }); + } + +private: + //! the class handling the primary variable switch + PriVarSwitchPtrTuple priVarSwitches_; + //! if we switched primary variables in the last iteration + std::array<bool, Assembler::Traits::numSubDomains> switchedInLastIteration_; +}; + +} // end namespace Dumux + +#endif