Line | Branch | Exec | Source |
---|---|---|---|
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-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
6 | // | ||
7 | /*! | ||
8 | * \file | ||
9 | * \ingroup OneEqModel | ||
10 | * | ||
11 | * \copydoc Dumux::OneEqVolumeVariables | ||
12 | */ | ||
13 | #ifndef DUMUX_ONEEQ_VOLUME_VARIABLES_HH | ||
14 | #define DUMUX_ONEEQ_VOLUME_VARIABLES_HH | ||
15 | |||
16 | #include <dune/common/math.hh> | ||
17 | #include <dumux/common/parameters.hh> | ||
18 | #include <dumux/freeflow/rans/volumevariables.hh> | ||
19 | |||
20 | namespace Dumux { | ||
21 | |||
22 | /*! | ||
23 | * \ingroup OneEqModel | ||
24 | * \brief Volume variables for the isothermal single-phase one-equation | ||
25 | * turbulence model by Spalart-Allmaras | ||
26 | */ | ||
27 | template <class Traits, class NSVolumeVariables> | ||
28 |
3/20✓ Branch 3 taken 74676 times.
✗ Branch 4 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 33032 times.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 33032 times.
✗ Branch 12 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
|
282360 | class OneEqVolumeVariables |
29 | : public RANSVolumeVariables<Traits, NSVolumeVariables> | ||
30 | { | ||
31 | using RANSParentType = RANSVolumeVariables<Traits, NSVolumeVariables>; | ||
32 | |||
33 | using Scalar = typename Traits::PrimaryVariables::value_type; | ||
34 | using DimVector = Dune::FieldVector<Scalar, Traits::ModelTraits::dim()>; | ||
35 | |||
36 | public: | ||
37 | //! export the indices type | ||
38 | using Indices = typename Traits::ModelTraits::Indices; | ||
39 | |||
40 | /*! | ||
41 | * \brief Update all quantities for a given control volume | ||
42 | * | ||
43 | * \param elemSol A vector containing all primary variables connected to the element | ||
44 | * \param problem The object specifying the problem which ought to | ||
45 | * be simulated | ||
46 | * \param element An element which contains part of the control volume | ||
47 | * \param scv The sub-control volume | ||
48 | */ | ||
49 | template<class ElementSolution, class Problem, class Element, class SubControlVolume> | ||
50 | 10449652 | void update(const ElementSolution &elemSol, | |
51 | const Problem &problem, | ||
52 | const Element &element, | ||
53 | const SubControlVolume& scv) | ||
54 | { | ||
55 | 10449652 | RANSParentType::updateNavierStokesVolVars(elemSol, problem, element, scv); | |
56 | 10449652 | updateRANSProperties(elemSol, problem, element, scv); | |
57 | 10449652 | } | |
58 | |||
59 | /*! | ||
60 | * \brief Update all turbulent quantities for a given control volume | ||
61 | * | ||
62 | * Wall and roughness related quantities are stored. Eddy viscosity is set. | ||
63 | * | ||
64 | * \param elemSol A vector containing all primary variables connected to the element | ||
65 | * \param problem The object specifying the problem which ought to be simulated | ||
66 | * \param element An element which contains part of the control volume | ||
67 | * \param scv The sub-control volume | ||
68 | */ | ||
69 | template<class ElementSolution, class Problem, class Element, class SubControlVolume> | ||
70 | 10449652 | void updateRANSProperties(const ElementSolution &elemSol, | |
71 | const Problem &problem, | ||
72 | const Element &element, | ||
73 | const SubControlVolume& scv) | ||
74 | { | ||
75 | 10449652 | RANSParentType::updateRANSProperties(elemSol, problem, element, scv); | |
76 | 10449652 | viscosityTilde_ = elemSol[0][Indices::viscosityTildeIdx]; | |
77 | 10449652 | storedViscosityTilde_ = problem.storedViscosityTilde(RANSParentType::elementIdx()); | |
78 | 10449652 | storedViscosityTildeGradient_ = problem.storedViscosityTildeGradient(RANSParentType::elementIdx()); | |
79 | 10449652 | stressTensorScalarProduct_ = problem.stressTensorScalarProduct(RANSParentType::elementIdx()); | |
80 | 10449652 | vorticityTensorScalarProduct_ = problem.vorticityTensorScalarProduct(RANSParentType::elementIdx()); | |
81 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 10449652 times.
|
10449652 | if (problem.useStoredEddyViscosity()) |
82 | ✗ | RANSParentType::setDynamicEddyViscosity_(problem.storedDynamicEddyViscosity(RANSParentType::elementIdx())); | |
83 | else | ||
84 | 20899304 | RANSParentType::setDynamicEddyViscosity_(calculateEddyViscosity()); | |
85 | 10449652 | RANSParentType::calculateEddyDiffusivity(problem); | |
86 | 10449652 | RANSParentType::calculateEddyThermalConductivity(problem); | |
87 | 10449652 | } | |
88 | |||
89 | /*! | ||
90 | * \brief Returns the dynamic eddy viscosity \f$\mathrm{[Pa s]}\f$. | ||
91 | */ | ||
92 | Scalar calculateEddyViscosity() | ||
93 |
1/8✓ Branch 0 taken 33032 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
10482684 | { return viscosityTilde() * fv1() * RANSParentType::density(); } |
94 | |||
95 | /*! | ||
96 | * \brief Returns the viscosity parameter \f$ m^2/s \f$ | ||
97 | */ | ||
98 | ✗ | Scalar viscosityTilde() const | |
99 | ✗ | { return viscosityTilde_; } | |
100 | |||
101 | /*! | ||
102 | * \brief Returns the viscosity parameter from the last iteration \f$ m^2/s \f$ | ||
103 | */ | ||
104 | Scalar storedViscosityTilde() const | ||
105 | { return storedViscosityTilde_; } | ||
106 | |||
107 | /*! | ||
108 | * \brief Returns the gradient of the viscosity parameter | ||
109 | */ | ||
110 | ✗ | DimVector storedViscosityTildeGradient() const | |
111 | ✗ | { return storedViscosityTildeGradient_; } | |
112 | |||
113 | /*! | ||
114 | * \brief Returns the scalar product of the stress tensor | ||
115 | */ | ||
116 | Scalar stressTensorScalarProduct() const | ||
117 | { return stressTensorScalarProduct_; } | ||
118 | |||
119 | /*! | ||
120 | * \brief Returns damping function for the eddy viscosity | ||
121 | */ | ||
122 | Scalar fv1() const | ||
123 | { | ||
124 | 153176528 | return viscosityRatio() * viscosityRatio() * viscosityRatio() | |
125 | 153176528 | / (viscosityRatio() * viscosityRatio() * viscosityRatio() | |
126 | 38294132 | + cv1() * cv1() * cv1()); | |
127 | } | ||
128 | |||
129 | //! \brief Returns a model function | ||
130 | Scalar fv2() const | ||
131 | 111245792 | { return 1.0 - viscosityRatio() / (1.0 + viscosityRatio() * fv1()); } | |
132 | |||
133 | //! \brief Returns a model function | ||
134 | ✗ | Scalar ft2() const | |
135 | { | ||
136 | using std::exp; | ||
137 | // the trip correction term is dropped according to Versteeg2009 and Wilcox2006 | ||
138 | // return ct3() * exp(-ct4() * viscosityRatio() * viscosityRatio()); | ||
139 | ✗ | return 0.0; | |
140 | } | ||
141 | |||
142 | //! \brief Returns a model function | ||
143 | 3973064 | Scalar fW() const | |
144 | { | ||
145 | using std::pow; | ||
146 | using Dune::power; | ||
147 | 7946128 | return g() * pow(1.0 + power(cw3(), 6) / (power(g(), 6) + power(cw3(), 6)), 1.0/6.0); | |
148 | } | ||
149 | |||
150 | //! \brief Returns a model function | ||
151 | 7946128 | Scalar g() const | |
152 | { | ||
153 | using Dune::power; | ||
154 | 15892256 | return r() + cw2() * (power(r(), 6) - r()); | |
155 | } | ||
156 | |||
157 | //! \brief Returns a model function | ||
158 | 23838384 | Scalar r() const | |
159 | { | ||
160 | using std::min; | ||
161 | 71515152 | return min(10.0, | |
162 | 23838384 | viscosityTilde() / stressTensorScalarProductTilde() | |
163 | 23838384 | / RANSParentType::karmanConstant() / RANSParentType::karmanConstant() | |
164 |
2/2✓ Branch 0 taken 23821428 times.
✓ Branch 1 taken 16956 times.
|
23838384 | / RANSParentType::wallDistance() / RANSParentType::wallDistance()); |
165 | } | ||
166 | |||
167 | //! \brief Returns the ratio of the kinematic viscosity and the viscosity parameter | ||
168 | Scalar viscosityRatio() const | ||
169 |
21/56✓ Branch 0 taken 150822 times.
✓ Branch 1 taken 27660626 times.
✓ Branch 2 taken 150822 times.
✓ Branch 3 taken 27660626 times.
✓ Branch 4 taken 150822 times.
✓ Branch 5 taken 27660626 times.
✓ Branch 6 taken 150822 times.
✓ Branch 7 taken 27660626 times.
✓ Branch 8 taken 150822 times.
✓ Branch 9 taken 27660626 times.
✓ Branch 10 taken 150822 times.
✓ Branch 11 taken 27660626 times.
✓ Branch 12 taken 150822 times.
✓ Branch 13 taken 27660626 times.
✓ Branch 14 taken 33032 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 33032 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 33032 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 33032 times.
✗ Branch 21 not taken.
✓ Branch 22 taken 33032 times.
✗ Branch 23 not taken.
✓ Branch 24 taken 33032 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 33032 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✗ Branch 47 not taken.
✗ Branch 48 not taken.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✗ Branch 51 not taken.
✗ Branch 52 not taken.
✗ Branch 53 not taken.
✗ Branch 54 not taken.
✗ Branch 55 not taken.
|
268058924 | { return viscosityTilde() / RANSParentType::kinematicViscosity(); } |
170 | |||
171 | /* | ||
172 | * ! \brief Returns a modified version of the stress tensor scalar product | ||
173 | * | ||
174 | * According to <a href="https://turbmodels.larc.nasa.gov/spalart.html">NASA</a> | ||
175 | * this term should never be zero and different limiters might be used. | ||
176 | * The Implementation uses the one proposed in: | ||
177 | * Allmaras, S. R., Johnson, F. T., and Spalart, P. R., | ||
178 | * "Modifications and Clarifications for the Implementation of the Spalart-Allmaras Turbulence Model," ICCFD7-1902 | ||
179 | */ | ||
180 | 27811448 | Scalar stressTensorScalarProductTilde() const | |
181 | { | ||
182 | // original form | ||
183 | // return vorticityMagnitude() | ||
184 | // + viscosityTilde() * fv2() | ||
185 | // / RANSParentType::karmanConstant() / RANSParentType::karmanConstant() | ||
186 | // / RANSParentType::wallDistance() / RANSParentType::wallDistance(); | ||
187 | |||
188 | // limiter form, literature source see above | ||
189 |
2/2✓ Branch 0 taken 150822 times.
✓ Branch 1 taken 27660626 times.
|
27811448 | Scalar sBar = viscosityTilde() * fv2() |
190 | 27811448 | / RANSParentType::karmanConstant() / RANSParentType::karmanConstant() | |
191 | 27811448 | / RANSParentType::wallDistance() / RANSParentType::wallDistance(); | |
192 | 55622896 | return sBar < -c2() * vorticityMagnitude() | |
193 | 55622896 | ? vorticityMagnitude() | |
194 | 603288 | + (vorticityMagnitude() * (c2() * c2() * vorticityMagnitude() + c3() * sBar)) | |
195 | 301644 | / ((c3() - 2.0 * c2()) * vorticityMagnitude() - sBar) | |
196 | 55472074 | : vorticityMagnitude() + sBar; | |
197 | } | ||
198 | |||
199 | //! \brief Returns the magnitude of the vorticity | ||
200 | ✗ | Scalar vorticityMagnitude() const | |
201 | { | ||
202 | using std::sqrt; | ||
203 |
2/2✓ Branch 0 taken 150822 times.
✓ Branch 1 taken 27660626 times.
|
27811448 | return sqrt(2.0 * vorticityTensorScalarProduct_); |
204 | } | ||
205 | |||
206 | //! \brief Returns a model constant | ||
207 | ✗ | Scalar c2() const | |
208 | ✗ | { return 0.7; } | |
209 | |||
210 | //! \brief Returns a model constant | ||
211 | ✗ | Scalar c3() const | |
212 | ✗ | { return 0.9; } | |
213 | |||
214 | //! \brief Returns a model constant | ||
215 | ✗ | Scalar sigma() const | |
216 | ✗ | { return 2.0/3.0; } | |
217 | |||
218 | //! \brief Returns a model constant | ||
219 | ✗ | Scalar cb1() const | |
220 | ✗ | { return 0.1355; } | |
221 | |||
222 | //! \brief Returns a model constant | ||
223 | ✗ | Scalar cb2() const | |
224 | ✗ | { return 0.622; } | |
225 | |||
226 | //! \brief Returns a model constant | ||
227 | ✗ | Scalar cv1() const | |
228 | ✗ | { return 7.1; } | |
229 | |||
230 | //! \brief Returns a model constant | ||
231 | Scalar ct3() const | ||
232 | { return 1.2; } | ||
233 | |||
234 | //! \brief Returns a model constant | ||
235 | Scalar ct4() const | ||
236 | { return 0.5; } | ||
237 | |||
238 | //! \brief Returns a model constant | ||
239 | Scalar cw1() const | ||
240 | { | ||
241 | 3973064 | return cb1() / RANSParentType::karmanConstant() / RANSParentType::karmanConstant() | |
242 | 3973064 | + (1.0 + cb2()) / sigma(); | |
243 | } | ||
244 | |||
245 | //! \brief Returns a model constant | ||
246 | ✗ | Scalar cw2() const | |
247 | ✗ | { return 0.3; } | |
248 | |||
249 | //! \brief Returns a model constant | ||
250 | ✗ | Scalar cw3() const | |
251 | ✗ | { return 2.0; } | |
252 | |||
253 | protected: | ||
254 | Scalar viscosityTilde_ = 0.0; | ||
255 | Scalar storedViscosityTilde_ = 0.0; | ||
256 | DimVector storedViscosityTildeGradient_ = DimVector(0.0); | ||
257 | Scalar stressTensorScalarProduct_; | ||
258 | Scalar vorticityTensorScalarProduct_; | ||
259 | }; | ||
260 | |||
261 | } // end namespace Dumux | ||
262 | |||
263 | #endif | ||
264 |