3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
staggered/elementsolution.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 *****************************************************************************/
24#ifndef DUMUX_STAGGERED_ELEMENT_SOLUTION_HH
25#define DUMUX_STAGGERED_ELEMENT_SOLUTION_HH
26
27#include <type_traits>
28#include <dune/istl/bvector.hh>
30
31namespace Dumux {
32
40template<class PrimaryVariables, class CellCenterPrimaryVariables>
41PrimaryVariables makePriVarsFromCellCenterPriVars(const CellCenterPrimaryVariables& cellCenterPriVars)
42{
43 static_assert(int(PrimaryVariables::dimension) > int(CellCenterPrimaryVariables::dimension),
44 "PrimaryVariables' size must be greater than the one of CellCenterPrimaryVariables");
45
46 PrimaryVariables priVars(0.0);
47 constexpr auto offset = PrimaryVariables::dimension - CellCenterPrimaryVariables::dimension;
48 for (std::size_t i = 0; i < cellCenterPriVars.size(); ++i)
49 priVars[i + offset] = cellCenterPriVars[i];
50 return priVars;
51}
52
53template<class PrimaryVariables>
54using StaggeredElementSolution = Dune::BlockVector<PrimaryVariables>;
55
61template<class FVElementGeometry, class PrimaryVariables>
62auto elementSolution(PrimaryVariables&& priVars)
63-> std::enable_if_t<FVElementGeometry::GridGeometry::discMethod == DiscretizationMethod::staggered,
65{
66 return StaggeredElementSolution<PrimaryVariables>({std::move(priVars)});
67}
68
76template<class PrimaryVariables, class CellCenterPrimaryVariables>
78{
79 return StaggeredElementSolution<PrimaryVariables>({makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars)});
80}
81
82} // end namespace Dumux
83
84#endif
The available discretization methods in Dumux.
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethod::box, BoxElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for box schemes.
Definition: box/elementsolution.hh:115
PrimaryVariables makePriVarsFromCellCenterPriVars(const CellCenterPrimaryVariables &cellCenterPriVars)
Helper function to create a PrimaryVariables object from CellCenterPrimaryVariables.
Definition: staggered/elementsolution.hh:41
StaggeredElementSolution< PrimaryVariables > makeElementSolutionFromCellCenterPrivars(const CellCenterPrimaryVariables &cellCenterPriVars)
Helper function to create an elementSolution from cell center primary variables.
Definition: staggered/elementsolution.hh:77
Definition: adapt.hh:29
Dune::BlockVector< PrimaryVariables > StaggeredElementSolution
Definition: staggered/elementsolution.hh:54