3.2-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
multidomain/newtonsolver.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
25#ifndef DUMUX_MULTIDOMAIN_NEWTON_SOLVER_HH
26#define DUMUX_MULTIDOMAIN_NEWTON_SOLVER_HH
27
28#include <memory>
30
31namespace Dumux {
32namespace Detail {
33
34template<class Assembler, class Index>
35using DetectPVSwitchMultiDomain = typename Assembler::template GridVariables<Index::value>::VolumeVariables::PrimaryVariableSwitch;
36
37template<class Assembler, std::size_t i>
38using GetPVSwitchMultiDomain = Dune::Std::detected_or<int, DetectPVSwitchMultiDomain, Assembler, Dune::index_constant<i>>;
39
40} // end namespace Detail
41
47template <class Assembler, class LinearSolver, class CouplingManager,
48 class Reassembler = DefaultPartialReassembler,
49 class Comm = Dune::CollectiveCommunication<Dune::MPIHelper::MPICommunicator> >
50class MultiDomainNewtonSolver: public NewtonSolver<Assembler, LinearSolver, Reassembler, Comm>
51{
53 using SolutionVector = typename Assembler::ResidualType;
54
55 template<std::size_t i>
56 using PrimaryVariableSwitch = typename Detail::GetPVSwitchMultiDomain<Assembler, i>::type;
57
58 template<std::size_t i>
59 using HasPriVarsSwitch = typename Detail::GetPVSwitchMultiDomain<Assembler, i>::value_t; // std::true_type or std::false_type
60
61 template<std::size_t i>
62 using PrivarSwitchPtr = std::unique_ptr<PrimaryVariableSwitch<i>>;
63 using PriVarSwitchPtrTuple = typename Assembler::Traits::template Tuple<PrivarSwitchPtr>;
64
65public:
66
70 MultiDomainNewtonSolver(std::shared_ptr<Assembler> assembler,
71 std::shared_ptr<LinearSolver> linearSolver,
72 std::shared_ptr<CouplingManager> couplingManager,
73 const Comm& comm = Dune::MPIHelper::getCollectiveCommunication(),
74 const std::string& paramGroup = "")
76 , couplingManager_(couplingManager)
77 {
78 using namespace Dune::Hybrid;
79 forEach(std::make_index_sequence<Assembler::Traits::numSubDomains>{}, [&](auto&& id)
80 {
81 const int priVarSwitchVerbosity = getParamFromGroup<int>(paramGroup, "PrimaryVariableSwitch.Verbosity", 1);
82 using PVSwitch = PrimaryVariableSwitch<std::decay_t<decltype(id)>::value>;
83 elementAt(priVarSwitches_, id) = std::make_unique<PVSwitch>(priVarSwitchVerbosity);
84 });
85
86 priVarsSwitchedInLastIteration_.fill(false);
87 }
88
92 void newtonBeginStep(const SolutionVector& uCurrentIter) override
93 {
94 ParentType::newtonBeginStep(uCurrentIter);
95 couplingManager_->updateSolution(uCurrentIter);
96 }
97
105 void newtonBegin(SolutionVector& u) override
106 {
108
109 using namespace Dune::Hybrid;
110 forEach(std::make_index_sequence<Assembler::Traits::numSubDomains>{}, [&](auto&& id)
111 {
112 this->initPriVarSwitch_(u, id, HasPriVarsSwitch<std::decay_t<decltype(id)>::value>{});
113 });
114 }
115
120 bool newtonConverged() const override
121 {
122 if (Dune::any_true(priVarsSwitchedInLastIteration_))
123 return false;
124
126 }
127
128
135 void newtonEndStep(SolutionVector &uCurrentIter, const SolutionVector &uLastIter) override
136 {
137 using namespace Dune::Hybrid;
138 forEach(std::make_index_sequence<Assembler::Traits::numSubDomains>{}, [&](auto&& id)
139 {
140 this->invokePriVarSwitch_(uCurrentIter[id], id, HasPriVarsSwitch<std::decay_t<decltype(id)>::value>{});
141 });
142
143 ParentType::newtonEndStep(uCurrentIter, uLastIter);
144 couplingManager_->updateSolution(uCurrentIter);
145 }
146
147private:
148
152 template<std::size_t i>
153 void initPriVarSwitch_(SolutionVector&, Dune::index_constant<i> id, std::false_type) {}
154
158 template<std::size_t i>
159 void initPriVarSwitch_(SolutionVector& sol, Dune::index_constant<i> id, std::true_type)
160 {
161 using namespace Dune::Hybrid;
162 auto& priVarSwitch = *elementAt(priVarSwitches_, id);
163
164 priVarSwitch.reset(sol[id].size());
165 priVarsSwitchedInLastIteration_[i] = false;
166
167 const auto& problem = this->assembler().problem(id);
168 const auto& gridGeometry = this->assembler().gridGeometry(id);
169 auto& gridVariables = this->assembler().gridVariables(id);
170 priVarSwitch.updateBoundary(problem, gridGeometry, gridVariables, sol[id]);
171 }
172
176 template<class SubSol, std::size_t i>
177 void invokePriVarSwitch_(SubSol&, Dune::index_constant<i> id, std::false_type) {}
178
182 template<class SubSol, std::size_t i>
183 void invokePriVarSwitch_(SubSol& uCurrentIter, Dune::index_constant<i> id, std::true_type)
184 {
185 // update the variable switch (returns true if the pri vars at at least one dof were switched)
186 // for disabled grid variable caching
187 const auto& gridGeometry = this->assembler().gridGeometry(id);
188 const auto& problem = this->assembler().problem(id);
189 auto& gridVariables = this->assembler().gridVariables(id);
190
191 using namespace Dune::Hybrid;
192 auto& priVarSwitch = *elementAt(priVarSwitches_, id);
193
194 // invoke the primary variable switch
195 priVarsSwitchedInLastIteration_[i] = priVarSwitch.update(uCurrentIter, gridVariables,
196 problem, gridGeometry);
197
198 if (priVarsSwitchedInLastIteration_[i])
199 {
200 for (const auto& element : elements(gridGeometry.gridView()))
201 {
202 // if the volume variables are cached globally, we need to update those where the primary variables have been switched
203 priVarSwitch.updateSwitchedVolVars(problem, element, gridGeometry, gridVariables, uCurrentIter);
204
205 // if the flux variables are cached globally, we need to update those where the primary variables have been switched
206 priVarSwitch.updateSwitchedFluxVarsCache(problem, element, gridGeometry, gridVariables, uCurrentIter);
207 }
208 }
209 }
210
212 std::shared_ptr<CouplingManager> couplingManager_;
213
215 PriVarSwitchPtrTuple priVarSwitches_;
217 std::array<bool, Assembler::Traits::numSubDomains> priVarsSwitchedInLastIteration_;
218};
219
220} // end namespace Dumux
221
222#endif
Definition: adapt.hh:29
typename Assembler::template GridVariables< Index::value >::VolumeVariables::PrimaryVariableSwitch DetectPVSwitchMultiDomain
Definition: multidomain/newtonsolver.hh:35
Dune::Std::detected_or< int, DetectPVSwitchMultiDomain, Assembler, Dune::index_constant< i > > GetPVSwitchMultiDomain
Definition: multidomain/newtonsolver.hh:38
Definition: partialreassembler.hh:42
const Assembler & assembler() const
Access the assembler.
Definition: common/pdesolver.hh:92
const LinearSolver & linearSolver() const
Access the linear solver.
Definition: common/pdesolver.hh:104
Base class for linear solvers.
Definition: solver.hh:37
Definition: multidomain/couplingmanager.hh:46
Newton sover for coupled problems.
Definition: multidomain/newtonsolver.hh:51
bool newtonConverged() const override
Returns true if the error of the solution is below the tolerance.
Definition: multidomain/newtonsolver.hh:120
void newtonEndStep(SolutionVector &uCurrentIter, const SolutionVector &uLastIter) override
Indicates that one Newton iteration was finished.
Definition: multidomain/newtonsolver.hh:135
void newtonBeginStep(const SolutionVector &uCurrentIter) override
Indicates the beginning of a Newton iteration.
Definition: multidomain/newtonsolver.hh:92
void newtonBegin(SolutionVector &u) override
Called before the Newton method is applied to an non-linear system of equations.
Definition: multidomain/newtonsolver.hh:105
MultiDomainNewtonSolver(std::shared_ptr< Assembler > assembler, std::shared_ptr< LinearSolver > linearSolver, std::shared_ptr< CouplingManager > couplingManager, const Comm &comm=Dune::MPIHelper::getCollectiveCommunication(), const std::string &paramGroup="")
The constructor.
Definition: multidomain/newtonsolver.hh:70
An implementation of a Newton solver.
Definition: nonlinear/newtonsolver.hh:97
const std::string & paramGroup() const
Returns the parameter group.
Definition: nonlinear/newtonsolver.hh:697
virtual void newtonBeginStep(const SolutionVector &u)
Indicates the beginning of a Newton iteration.
Definition: nonlinear/newtonsolver.hh:320
virtual void newtonEndStep(SolutionVector &uCurrentIter, const SolutionVector &uLastIter)
Indicates that one Newton iteration was finished.
Definition: nonlinear/newtonsolver.hh:493
const Communication & comm() const
the communicator for parallel runs
Definition: nonlinear/newtonsolver.hh:146
virtual void newtonBegin(SolutionVector &u)
Called before the Newton method is applied to an non-linear system of equations.
Definition: nonlinear/newtonsolver.hh:270
virtual bool newtonConverged() const
Returns true if the error of the solution is below the tolerance.
Definition: nonlinear/newtonsolver.hh:541
Reference implementation of the Newton solver.