3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
maxwellstefandiffusioncoefficients.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 *****************************************************************************/
26#ifndef DUMUX_FLUX_MAXWELLSTEFAN_DIFFUSION_COEFFICIENTS_HH
27#define DUMUX_FLUX_MAXWELLSTEFAN_DIFFUSION_COEFFICIENTS_HH
28
29#include <array>
30#include <cassert>
31
32namespace Dumux {
33
43template <class Scalar, int numPhases, int numComponents>
45{
46public:
47 template<class DiffCoeffFunc>
48 void update(const DiffCoeffFunc& computeDiffCoeff)
49 {
50 for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx)
51 for (int compIIdx = 0; compIIdx < numComponents; ++compIIdx)
52 for (int compJIdx = compIIdx+1; compJIdx < numComponents; ++compJIdx)
53 diffCoeff_[getIndex_(phaseIdx, compIIdx, compJIdx)]
54 = computeDiffCoeff(phaseIdx, compIIdx, compJIdx);
55 }
56
57 Scalar operator() (int phaseIdx, int compIIdx, int compJIdx) const
58 {
59 sortComponentIndices_(compIIdx, compJIdx);
60 assert(compIIdx != compJIdx);
61 return diffCoeff_[getIndex_(phaseIdx, compIIdx, compJIdx)];
62 }
63
64private:
75 std::array<Scalar, (numPhases * ((numComponents * (numComponents - 1)) / 2))> diffCoeff_;
76
88 constexpr int getIndex_(int phaseIdx, int compIIdx, int compJIdx) const
89 {
90 return phaseIdx * ((numComponents * (numComponents - 1)) / 2)
91 + compIIdx * numComponents
92 - ((compIIdx * (compIIdx + 1)) / 2)
93 + compJIdx - (compIIdx +1);
94 }
95
96 void sortComponentIndices_(int& compIIdx, int& compJIdx) const
97 { if (compIIdx > compJIdx) std::swap(compIIdx, compJIdx); }
98};
99
100} // end namespace Dumux
101
102#endif
Definition: adapt.hh:29
Container storing the diffusion coefficients required by the Maxwell- Stefan diffusion law....
Definition: maxwellstefandiffusioncoefficients.hh:45
void update(const DiffCoeffFunc &computeDiffCoeff)
Definition: maxwellstefandiffusioncoefficients.hh:48
Scalar operator()(int phaseIdx, int compIIdx, int compJIdx) const
Definition: maxwellstefandiffusioncoefficients.hh:57