version 3.10.0
Loading...
Searching...
No Matches
multidomain/fvassembler.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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
14#ifndef DUMUX_MULTIDOMAIN_FV_ASSEMBLER_HH
15#define DUMUX_MULTIDOMAIN_FV_ASSEMBLER_HH
16
17#include <type_traits>
18#include <tuple>
19
20#include <dune/common/hybridutilities.hh>
21#include <dune/istl/matrixindexset.hh>
22
34
40#include "assemblerview.hh"
41
43
44namespace Dumux {
45
46namespace Grid::Capabilities {
47
48namespace Detail {
49// helper for multi-domain models
50template<class T, std::size_t... I>
51bool allGridsSupportsMultithreadingImpl(const T& gridGeometries, std::index_sequence<I...>)
52{
53 return (... && supportsMultithreading(std::get<I>(gridGeometries)->gridView()));
54}
55} // end namespace Detail
56
57// helper for multi-domain models (all grids have to support multithreading)
58template<class... GG>
59bool allGridsSupportsMultithreading(const std::tuple<GG...>& gridGeometries)
60{
61 return Detail::allGridsSupportsMultithreadingImpl<std::tuple<GG...>>(gridGeometries, std::make_index_sequence<sizeof...(GG)>());
62}
63
64} // end namespace Grid::Capabilities
65
71template<class CM>
72struct CouplingManagerSupportsMultithreadedAssembly : public std::false_type
73{};
74
84template<class MDTraits, class CMType, DiffMethod diffMethod, bool useImplicitAssembly = true>
86{
87 template<std::size_t id>
88 using SubDomainTypeTag = typename MDTraits::template SubDomain<id>::TypeTag;
89
90public:
91 using Traits = MDTraits;
92
93 using Scalar = typename MDTraits::Scalar;
94
95 template<std::size_t id>
96 using LocalResidual = typename MDTraits::template SubDomain<id>::LocalResidual;
97
98 template<std::size_t id>
99 using GridVariables = typename MDTraits::template SubDomain<id>::GridVariables;
100
101 template<std::size_t id>
102 using GridGeometry = typename MDTraits::template SubDomain<id>::GridGeometry;
103
104 template<std::size_t id>
105 using Problem = typename MDTraits::template SubDomain<id>::Problem;
106
107 using JacobianMatrix = typename MDTraits::JacobianMatrix;
108 using SolutionVector = typename MDTraits::SolutionVector;
109 using ResidualType = typename MDTraits::ResidualVector;
110
111 using CouplingManager = CMType;
112
116 static constexpr bool isImplicit()
117 { return useImplicitAssembly; }
118
119private:
120
121 using ProblemTuple = typename MDTraits::template TupleOfSharedPtrConst<Problem>;
122 using GridGeometryTuple = typename MDTraits::template TupleOfSharedPtrConst<GridGeometry>;
123 using GridVariablesTuple = typename MDTraits::template TupleOfSharedPtr<GridVariables>;
124
126 using ThisType = MultiDomainFVAssembler<MDTraits, CouplingManager, diffMethod, isImplicit()>;
127
128 template<std::size_t id>
129 using SubDomainAssemblerView = MultiDomainAssemblerSubDomainView<ThisType, id>;
130
131 template<class DiscretizationMethod, std::size_t id>
132 struct SubDomainAssemblerType;
133
134 template<std::size_t id>
135 struct SubDomainAssemblerType<DiscretizationMethods::CCTpfa, id>
136 {
137 using type = SubDomainCCLocalAssembler<id, SubDomainTypeTag<id>, SubDomainAssemblerView<id>, diffMethod, isImplicit()>;
138 };
139
140 template<std::size_t id>
141 struct SubDomainAssemblerType<DiscretizationMethods::CCMpfa, id>
142 {
143 using type = SubDomainCCLocalAssembler<id, SubDomainTypeTag<id>, SubDomainAssemblerView<id>, diffMethod, isImplicit()>;
144 };
145
146 template<std::size_t id, class DM>
147 struct SubDomainAssemblerType<DiscretizationMethods::CVFE<DM>, id>
148 {
149 using type = SubDomainCVFELocalAssembler<id, SubDomainTypeTag<id>, SubDomainAssemblerView<id>, diffMethod, isImplicit()>;
150 };
151
152 template<std::size_t id>
153 struct SubDomainAssemblerType<DiscretizationMethods::Staggered, id>
154 {
155 using type = SubDomainStaggeredLocalAssembler<id, SubDomainTypeTag<id>, SubDomainAssemblerView<id>, diffMethod, isImplicit()>;
156 };
157
158 template<std::size_t id>
159 struct SubDomainAssemblerType<DiscretizationMethods::FCStaggered, id>
160 {
161 using type = SubDomainFaceCenteredLocalAssembler<id, SubDomainTypeTag<id>, SubDomainAssemblerView<id>, diffMethod, isImplicit()>;
162 };
163
164 template<std::size_t id>
165 using SubDomainAssembler = typename SubDomainAssemblerType<typename GridGeometry<id>::DiscretizationMethod, id>::type;
166
167public:
168
169
176 GridGeometryTuple gridGeometry,
177 GridVariablesTuple gridVariables,
178 std::shared_ptr<CouplingManager> couplingManager)
180 , problemTuple_(std::move(problem))
181 , gridGeometryTuple_(std::move(gridGeometry))
182 , gridVariablesTuple_(std::move(gridVariables))
183 , timeLoop_()
184 , isStationaryProblem_(true)
185 , warningIssued_(false)
186 {
187 static_assert(isImplicit(), "Explicit assembler for stationary problem doesn't make sense!");
188 std::cout << "Instantiated assembler for a stationary problem." << std::endl;
189
193 && getParam<bool>("Assembly.Multithreading", true);
194
195 maybeComputeColors_();
196 }
197
204 GridGeometryTuple gridGeometry,
205 GridVariablesTuple gridVariables,
206 std::shared_ptr<CouplingManager> couplingManager,
207 std::shared_ptr<const TimeLoop> timeLoop,
208 const SolutionVector& prevSol)
210 , problemTuple_(std::move(problem))
211 , gridGeometryTuple_(std::move(gridGeometry))
212 , gridVariablesTuple_(std::move(gridVariables))
213 , timeLoop_(timeLoop)
214 , prevSol_(&prevSol)
215 , isStationaryProblem_(false)
216 , warningIssued_(false)
217 {
218 std::cout << "Instantiated assembler for an instationary problem." << std::endl;
219
223 && getParam<bool>("Assembly.Multithreading", true);
224
225 maybeComputeColors_();
226 }
227
233 {
234 checkAssemblerState_();
235 resetJacobian_();
236 resetResidual_();
237
238 using namespace Dune::Hybrid;
239 forEach(std::make_index_sequence<JacobianMatrix::N()>(), [&](const auto domainId)
240 {
241 auto& jacRow = (*jacobian_)[domainId];
242 auto& subRes = (*residual_)[domainId];
243 this->assembleJacobianAndResidual_(domainId, jacRow, subRes, curSol);
244
245 const auto gridGeometry = std::get<domainId>(gridGeometryTuple_);
246 enforcePeriodicConstraints_(domainId, jacRow, subRes, *gridGeometry, curSol[domainId]);
247 });
248 }
249
252 {
253 resetResidual_();
254 assembleResidual(*residual_, curSol);
255 }
256
259 {
260 r = 0.0;
261
262 checkAssemblerState_();
263
264 // update the grid variables for the case of active caching
265 updateGridVariables(curSol);
266
267 using namespace Dune::Hybrid;
268 forEach(integralRange(Dune::Hybrid::size(r)), [&](const auto domainId)
269 {
270 auto& subRes = r[domainId];
271 this->assembleResidual_(domainId, subRes, curSol);
272 });
273 }
274
280 void setLinearSystem(std::shared_ptr<JacobianMatrix> A,
281 std::shared_ptr<ResidualType> r)
282 {
283 jacobian_ = A;
284 residual_ = r;
285
286 setJacobianBuildMode(*jacobian_);
287 setJacobianPattern_(*jacobian_);
288 setResidualSize_(*residual_);
289 }
290
296 {
297 jacobian_ = std::make_shared<JacobianMatrix>();
298 residual_ = std::make_shared<ResidualType>();
299
300 setJacobianBuildMode(*jacobian_);
301 setJacobianPattern_(*jacobian_);
302 setResidualSize_(*residual_);
303 }
304
309 {
310 using namespace Dune::Hybrid;
311 forEach(std::make_index_sequence<JacobianMatrix::N()>(), [&](const auto i)
312 {
313 forEach(jac[i], [&](auto& jacBlock)
314 {
315 using BlockType = std::decay_t<decltype(jacBlock)>;
316 if (jacBlock.buildMode() == BlockType::BuildMode::unknown)
317 jacBlock.setBuildMode(BlockType::BuildMode::random);
318 else if (jacBlock.buildMode() != BlockType::BuildMode::random)
319 DUNE_THROW(Dune::NotImplemented, "Only BCRS matrices with random build mode are supported at the moment");
320 });
321 });
322 }
323
328 {
329 setJacobianPattern_(*jacobian_);
330 setResidualSize_(*residual_);
331 maybeComputeColors_();
332 }
333
338 {
339 using namespace Dune::Hybrid;
340 forEach(integralRange(Dune::Hybrid::size(gridVariablesTuple_)), [&](const auto domainId)
341 { this->gridVariables(domainId).update(curSol[domainId]); });
342 }
343
347 void resetTimeStep(const SolutionVector& curSol)
348 {
349 using namespace Dune::Hybrid;
350 forEach(integralRange(Dune::Hybrid::size(gridVariablesTuple_)), [&](const auto domainId)
351 { this->gridVariables(domainId).resetTimeStep(curSol[domainId]); });
352 }
353
355 template<std::size_t i>
356 std::size_t numDofs(Dune::index_constant<i> domainId) const
357 { return std::get<domainId>(gridGeometryTuple_)->numDofs(); }
358
360 template<std::size_t i>
361 const auto& problem(Dune::index_constant<i> domainId) const
362 { return *std::get<domainId>(problemTuple_); }
363
365 template<std::size_t i>
366 const auto& gridGeometry(Dune::index_constant<i> domainId) const
367 { return *std::get<domainId>(gridGeometryTuple_); }
368
370 template<std::size_t i>
371 const auto& gridView(Dune::index_constant<i> domainId) const
372 { return gridGeometry(domainId).gridView(); }
373
375 template<std::size_t i>
376 GridVariables<i>& gridVariables(Dune::index_constant<i> domainId)
377 { return *std::get<domainId>(gridVariablesTuple_); }
378
380 template<std::size_t i>
381 const GridVariables<i>& gridVariables(Dune::index_constant<i> domainId) const
382 { return *std::get<domainId>(gridVariablesTuple_); }
383
386 { return *couplingManager_; }
387
390 { return *jacobian_; }
391
394 { return *residual_; }
395
397 const SolutionVector& prevSol() const
398 { return *prevSol_; }
399
404 void setTimeManager(std::shared_ptr<const TimeLoop> timeLoop)
405 { timeLoop_ = timeLoop; isStationaryProblem_ = !(static_cast<bool>(timeLoop)); }
406
412 { prevSol_ = &u; }
413
418 { return isStationaryProblem_; }
419
423 template<std::size_t i>
424 LocalResidual<i> localResidual(Dune::index_constant<i> domainId) const
425 { return LocalResidual<i>(std::get<domainId>(problemTuple_).get(), timeLoop_.get()); }
426
427protected:
429 std::shared_ptr<CouplingManager> couplingManager_;
430
431private:
435 void setJacobianPattern_(JacobianMatrix& jac) const
436 {
437 using namespace Dune::Hybrid;
438 forEach(std::make_index_sequence<JacobianMatrix::N()>(), [&](const auto domainI)
439 {
440 forEach(integralRange(Dune::Hybrid::size(jac[domainI])), [&](const auto domainJ)
441 {
442 const auto pattern = this->getJacobianPattern_(domainI, domainJ);
443 pattern.exportIdx(jac[domainI][domainJ]);
444 });
445 });
446 }
447
451 void setResidualSize_(ResidualType& res) const
452 {
453 using namespace Dune::Hybrid;
454 forEach(integralRange(Dune::Hybrid::size(res)), [&](const auto domainId)
455 { res[domainId].resize(this->numDofs(domainId)); });
456 }
457
458 // reset the residual vector to 0.0
459 void resetResidual_()
460 {
461 if(!residual_)
462 {
463 residual_ = std::make_shared<ResidualType>();
464 setResidualSize_(*residual_);
465 }
466
467 (*residual_) = 0.0;
468 }
469
470 // reset the jacobian vector to 0.0
471 void resetJacobian_()
472 {
473 if(!jacobian_)
474 {
475 jacobian_ = std::make_shared<JacobianMatrix>();
476 setJacobianBuildMode(*jacobian_);
477 setJacobianPattern_(*jacobian_);
478 }
479
480 (*jacobian_) = 0.0;
481 }
482
484 void maybeComputeColors_()
485 {
486 if constexpr (CouplingManagerSupportsMultithreadedAssembly<CouplingManager>::value)
487 if (enableMultithreading_)
488 couplingManager_->computeColorsForAssembly();
489 }
490
491 // check if the assembler is in a correct state for assembly
492 void checkAssemblerState_() const
493 {
494 if (!isStationaryProblem_ && !prevSol_)
495 DUNE_THROW(Dune::InvalidStateException, "Assembling instationary problem but previous solution was not set!");
496
497 if (isStationaryProblem_ && prevSol_)
498 DUNE_THROW(Dune::InvalidStateException, "Assembling stationary problem but a previous solution was set."
499 << " Did you forget to set the timeLoop to make this problem instationary?");
500 }
501
502 template<std::size_t i, class JacRow, class SubRes>
503 void assembleJacobianAndResidual_(Dune::index_constant<i> domainId, JacRow& jacRow, SubRes& subRes,
504 const SolutionVector& curSol)
505 {
506 assemble_(domainId, [&](const auto& element)
507 {
508 MultiDomainAssemblerSubDomainView view{*this, domainId};
509 SubDomainAssembler<i> subDomainAssembler(view, element, curSol, *couplingManager_);
510 subDomainAssembler.assembleJacobianAndResidual(jacRow, subRes, gridVariablesTuple_);
511 });
512 }
513
514 template<std::size_t i, class SubRes>
515 void assembleResidual_(Dune::index_constant<i> domainId, SubRes& subRes,
516 const SolutionVector& curSol)
517 {
518 assemble_(domainId, [&](const auto& element)
519 {
520 MultiDomainAssemblerSubDomainView view{*this, domainId};
521 SubDomainAssembler<i> subDomainAssembler(view, element, curSol, *couplingManager_);
522 subDomainAssembler.assembleResidual(subRes);
523 });
524 }
525
531 template<std::size_t i, class AssembleElementFunc>
532 void assemble_(Dune::index_constant<i> domainId, AssembleElementFunc&& assembleElement) const
533 {
534 // a state that will be checked on all processes
535 bool succeeded = false;
536
537 // try assembling using the local assembly function
538 try
539 {
540 if constexpr (CouplingManagerSupportsMultithreadedAssembly<CouplingManager>::value)
541 {
542 if (enableMultithreading_)
543 {
544 couplingManager_->assembleMultithreaded(
545 domainId, std::forward<AssembleElementFunc>(assembleElement)
546 );
547 return;
548 }
549 }
550
551 // fallback for coupling managers that don't support multithreaded assembly (yet)
552 // or if multithreaded assembly is disabled
553 // let the local assembler add the element contributions
554 for (const auto& element : elements(gridView(domainId)))
555 assembleElement(element);
556
557 // if we get here, everything worked well on this process
558 succeeded = true;
559 }
560 // throw exception if a problem occurred
561 catch (NumericalProblem &e)
562 {
563 std::cout << "rank " << gridView(domainId).comm().rank()
564 << " caught an exception while assembling:" << e.what()
565 << "\n";
566 succeeded = false;
567 }
568
569 // make sure everything worked well on all processes
570 if (gridView(domainId).comm().size() > 1)
571 succeeded = gridView(domainId).comm().min(succeeded);
572
573 // if not succeeded rethrow the error on all processes
574 if (!succeeded)
575 DUNE_THROW(NumericalProblem, "A process did not succeed in linearizing the system");
576 }
577
578 // get diagonal block pattern
579 template<std::size_t i, std::size_t j, typename std::enable_if_t<(i==j), int> = 0>
580 Dune::MatrixIndexSet getJacobianPattern_(Dune::index_constant<i> domainI,
581 Dune::index_constant<j> domainJ) const
582 {
583 const auto& gg = gridGeometry(domainI);
584 auto pattern = getJacobianPattern<isImplicit()>(gg);
585 couplingManager_->extendJacobianPattern(domainI, pattern);
586 return pattern;
587 }
588
589 // get coupling block pattern
590 template<std::size_t i, std::size_t j, typename std::enable_if_t<(i!=j), int> = 0>
591 Dune::MatrixIndexSet getJacobianPattern_(Dune::index_constant<i> domainI,
592 Dune::index_constant<j> domainJ) const
593 {
595 domainI, gridGeometry(domainI),
596 domainJ, gridGeometry(domainJ));
597 }
598
599 // build periodic constraints into the system matrix
600 template<std::size_t i, class JacRow, class Res, class GG, class Sol>
601 void enforcePeriodicConstraints_(Dune::index_constant<i> domainI, JacRow& jacRow, Res& res, const GG& gridGeometry, const Sol& curSol)
602 {
603 if constexpr (Detail::hasPeriodicDofMap<GG>())
604 {
605 for (const auto& m : gridGeometry.periodicDofMap())
606 {
607 if (m.first < m.second)
608 {
609 auto& jac = jacRow[domainI];
610
611 // add the second row to the first
612 res[m.first] += res[m.second];
613
614 const auto end = jac[m.second].end();
615 for (auto it = jac[m.second].begin(); it != end; ++it)
616 jac[m.first][it.index()] += (*it);
617
618
619 // enforce the solution of the first periodic DOF to the second one
620 res[m.second] = curSol[m.second] - curSol[m.first];
621
622 // set derivatives accordingly in jacobian, i.e. id for m.second and -id for m.first
623 auto setMatrixBlock = [] (auto& matrixBlock, double diagValue)
624 {
625 for (int eIdx = 0; eIdx < matrixBlock.N(); ++eIdx)
626 matrixBlock[eIdx][eIdx] = diagValue;
627 };
628
629 for (auto it = jac[m.second].begin(); it != end; ++it)
630 {
631 auto& matrixBlock = *it;
632 matrixBlock = 0.0;
633
634 assert(matrixBlock.N() == matrixBlock.M());
635 if(it.index() == m.second)
636 setMatrixBlock(matrixBlock, 1.0);
637
638 if(it.index() == m.first)
639 setMatrixBlock(matrixBlock, -1.0);
640
641 }
642
643 using namespace Dune::Hybrid;
644 forEach(makeIncompleteIntegerSequence<JacRow::size(), domainI>(), [&](const auto couplingDomainId)
645 {
646 auto& jacCoupling = jacRow[couplingDomainId];
647
648 for (auto it = jacCoupling[m.second].begin(); it != jacCoupling[m.second].end(); ++it)
649 jacCoupling[m.first][it.index()] += (*it);
650
651 for (auto it = jacCoupling[m.second].begin(); it != jacCoupling[m.second].end(); ++it)
652 (*it) = 0.0;
653 });
654 }
655 }
656 }
657 }
658
660 ProblemTuple problemTuple_;
661
663 GridGeometryTuple gridGeometryTuple_;
664
666 GridVariablesTuple gridVariablesTuple_;
667
669 std::shared_ptr<const TimeLoop> timeLoop_;
670
672 const SolutionVector* prevSol_ = nullptr;
673
675 bool isStationaryProblem_;
676
678 std::shared_ptr<JacobianMatrix> jacobian_;
679 std::shared_ptr<ResidualType> residual_;
680
682 mutable bool warningIssued_;
683
685 bool enableMultithreading_ = false;
686};
687
688} // end namespace Dumux
689
690#endif
Subdomain-specific views on multidomain assemblers.
The interface of the coupling manager for multi domain problems.
Definition multidomain/couplingmanager.hh:37
Subdomain-specific view on a multidomain assembler. Allows retrieval of sub-domain specific objects w...
Definition assemblerview.hh:31
A linear system assembler (residual and Jacobian) for finite volume schemes (box, tpfa,...
Definition multidomain/fvassembler.hh:86
std::size_t numDofs(Dune::index_constant< i > domainId) const
the number of dof locations of domain i
Definition multidomain/fvassembler.hh:356
void updateGridVariables(const SolutionVector &curSol)
Updates the grid variables with the given solution.
Definition multidomain/fvassembler.hh:337
void updateAfterGridAdaption()
Resizes jacobian and residual and recomputes colors.
Definition multidomain/fvassembler.hh:327
typename StaggeredMultiDomainTraits< TypeTag, TypeTag >::template SubDomain< id >::Problem Problem
Definition multidomain/fvassembler.hh:105
JacobianMatrix & jacobian()
the full Jacobian matrix
Definition multidomain/fvassembler.hh:389
static constexpr bool isImplicit()
Returns true if the assembler considers implicit assembly.
Definition multidomain/fvassembler.hh:116
void assembleResidual(ResidualType &r, const SolutionVector &curSol)
assemble a residual r
Definition multidomain/fvassembler.hh:258
StaggeredCouplingManager< StaggeredMultiDomainTraits< TypeTag, TypeTag > > CouplingManager
Definition multidomain/fvassembler.hh:111
void assembleResidual(const SolutionVector &curSol)
compute the residuals using the internal residual
Definition multidomain/fvassembler.hh:251
void resetTimeStep(const SolutionVector &curSol)
Resets the grid variables to the last time step.
Definition multidomain/fvassembler.hh:347
const GridVariables< i > & gridVariables(Dune::index_constant< i > domainId) const
the grid variables of domain i
Definition multidomain/fvassembler.hh:381
void setPreviousSolution(const SolutionVector &u)
Sets the solution from which to start the time integration. Has to be called prior to assembly for ti...
Definition multidomain/fvassembler.hh:411
typename StaggeredMultiDomainTraits< TypeTag, TypeTag >::template SubDomain< id >::GridVariables GridVariables
Definition multidomain/fvassembler.hh:99
void setLinearSystem()
The version without arguments uses the default constructor to create the jacobian and residual object...
Definition multidomain/fvassembler.hh:295
const auto & gridView(Dune::index_constant< i > domainId) const
the grid view of domain i
Definition multidomain/fvassembler.hh:371
typename StaggeredMultiDomainTraits< TypeTag, TypeTag >::template SubDomain< id >::GridGeometry GridGeometry
Definition multidomain/fvassembler.hh:102
void setTimeManager(std::shared_ptr< const TimeLoop > timeLoop)
Set time loop for instationary problems.
Definition multidomain/fvassembler.hh:404
MultiDomainFVAssembler(ProblemTuple problem, GridGeometryTuple gridGeometry, GridVariablesTuple gridVariables, std::shared_ptr< CouplingManager > couplingManager, std::shared_ptr< const TimeLoop > timeLoop, const SolutionVector &prevSol)
The constructor for instationary problems.
Definition multidomain/fvassembler.hh:203
void setJacobianBuildMode(JacobianMatrix &jac) const
Sets the jacobian build mode.
Definition multidomain/fvassembler.hh:308
typename StaggeredMultiDomainTraits< TypeTag, TypeTag >::template SubDomain< id >::LocalResidual LocalResidual
Definition multidomain/fvassembler.hh:96
LocalResidual< i > localResidual(Dune::index_constant< i > domainId) const
Create a local residual object (used by the local assembler)
Definition multidomain/fvassembler.hh:424
MultiDomainFVAssembler(ProblemTuple problem, GridGeometryTuple gridGeometry, GridVariablesTuple gridVariables, std::shared_ptr< CouplingManager > couplingManager)
The constructor for stationary problems.
Definition multidomain/fvassembler.hh:175
void assembleJacobianAndResidual(const SolutionVector &curSol)
Assembles the global Jacobian of the residual and the residual for the current solution.
Definition multidomain/fvassembler.hh:232
bool isStationaryProblem() const
Whether we are assembling a stationary or instationary problem.
Definition multidomain/fvassembler.hh:417
ResidualType & residual()
the full residual vector
Definition multidomain/fvassembler.hh:393
void setLinearSystem(std::shared_ptr< JacobianMatrix > A, std::shared_ptr< ResidualType > r)
Tells the assembler which jacobian and residual to use. This also resizes the containers to the requi...
Definition multidomain/fvassembler.hh:280
The cell-centered scheme multidomain local assembler.
Definition multidomain/subdomaincclocalassembler.hh:270
Manages the handling of time dependent problems.
Definition common/timeloop.hh:84
The default time loop for instationary simulations.
Definition common/timeloop.hh:139
Defines all properties used in Dumux.
Manages the handling of time dependent problems.
Helper function to generate Jacobian pattern for multi domain models.
An enum class to define various differentiation methods available in order to compute the derivatives...
Some exceptions thrown in DuMux
dune-grid capabilities compatibility layer
Dune::MatrixIndexSet getJacobianPattern(const GridGeometry &gridGeometry)
Helper function to generate Jacobian pattern for cell-centered methods.
Definition jacobianpattern.hh:28
Dune::MatrixIndexSet getCouplingJacobianPattern(const CouplingManager &couplingManager, Dune::index_constant< i > domainI, const GridGeometryI &gridGeometryI, Dune::index_constant< j > domainJ, const GridGeometryJ &gridGeometryJ)
Helper function to generate coupling Jacobian pattern (off-diagonal blocks) for cell-centered schemes...
Definition couplingjacobianpattern.hh:30
constexpr bool isSerial()
Checking whether the backend is serial.
Definition multithreading.hh:45
T getParam(Args &&... args)
A free function to get a parameter from the parameter tree singleton.
Definition parameters.hh:139
Helper function to generate Jacobian pattern for different discretization methods.
The available discretization methods in Dumux.
A multidomain local assembler for Jacobian and residual contribution per element (cell-centered metho...
An assembler for Jacobian and residual contribution per element (CVFE methods) for multidomain proble...
Multithreading in Dumux.
constexpr bool hasPeriodicDofMap()
Definition periodic.hh:24
Definition method.hh:20
Definition multistagemultidomainfvassembler.hh:49
bool allGridsSupportsMultithreadingImpl(const T &gridGeometries, std::index_sequence< I... >)
Definition multistagemultidomainfvassembler.hh:52
Definition gridcapabilities.hh:57
bool allGridsSupportsMultithreading(const std::tuple< GG... > &gridGeometries)
Definition multistagemultidomainfvassembler.hh:60
bool supportsMultithreading(const GridView &gridView)
Definition gridcapabilities.hh:74
Definition adapt.hh:17
const Scalar PengRobinsonMixture< Scalar, StaticParameters >::u
Definition pengrobinsonmixture.hh:138
typename Detail::ConcatSeq< decltype(std::make_index_sequence< e >{}), e+1, decltype(std::make_index_sequence<(n > e) ?(n - e - 1) :0 >{})>::type makeIncompleteIntegerSequence
Definition utility.hh:58
Provides a helper class for nonoverlapping decomposition.
Type traits to detect periodicity support.
Type trait that is specialized for coupling manager supporting multithreaded assembly.
Definition multistagemultidomainfvassembler.hh:78
An assembler for Jacobian and residual contribution per element (face-centered staggered methods) for...
A multidomain assembler for Jacobian and residual contribution per element (staggered method)
Utilities for template meta programming.