OpenLB 1.8.1
Loading...
Searching...
No Matches
collisionLES.h
Go to the documentation of this file.
1/* This file is part of the OpenLB library
2 *
3 * Copyright (C) 2012-2015 Mathias J. Krause, Jonas Latt, Patrick Nathen
4 * 2021 Adrian Kummerlaender
5 * E-mail contact: info@openlb.net
6 * The most recent release of OpenLB can be downloaded at
7 * <http://www.openlb.net/>
8 *
9 * This program is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU General Public License
11 * as published by the Free Software Foundation; either version 2
12 * of the License, or (at your option) any later version.
13 *
14 * This program is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 * GNU General Public License for more details.
18 *
19 * You should have received a copy of the GNU General Public
20 * License along with this program; if not, write to the Free
21 * Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
22 * Boston, MA 02110-1301, USA.
23*/
24
25#ifndef DYNAMICS_COLLISION_LES_H
26#define DYNAMICS_COLLISION_LES_H
27
28#include "collision.h"
29#include "collisionModifiers.h"
30
31namespace olb {
32
33namespace collision {
34
35namespace LES {
36
37struct SMAGORINSKY : public descriptors::FIELD_BASE<1,0,0> {
38 template <typename T, typename DESCRIPTOR>
39 static constexpr auto getInitialValue() {
40 return FieldD<T,DESCRIPTOR,SMAGORINSKY>(0.15);
41 }
42
43 template <typename T, typename DESCRIPTOR,typename FIELD>
44 static constexpr auto isValid(FieldD<T,DESCRIPTOR,FIELD> value) {
45 return value >= 0;
46 }
47};
48
49
50}
51
53
56namespace detail {
57
58template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
59struct SmagorinskyEffectiveOmega {
60 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
61 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
62
63 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
64 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
65 V piNeqNormSqr { };
66 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
67 const V rho = MomentaF().computeRho(cell);
68 const V omega = parameters.template get<descriptors::OMEGA>();
69 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
70 V piNeqNorm = util::sqrt(piNeqNormSqr);
71 V preFactor = smagorinsky*smagorinsky
72 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
73 * 2 * util::sqrt(2);
75 V tauMol = V{1} / omega;
77 V tauTurb = V{0.5} * (util::sqrt(tauMol*tauMol + preFactor / rho * piNeqNorm) - tauMol);
79 V tauEff = tauMol + tauTurb;
80 return V{1} / tauEff;
81 }
82
83 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
84 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
85 parameters.template set<descriptors::OMEGA>(
86 computeEffectiveOmega(cell, parameters));
87 return CollisionO().apply(cell, parameters);
88 }
89};
90
91template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
92struct IncompressibleSmagorinskyEffectiveOmega {
93 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
94 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
95
96 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
97 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
98 V piNeqNormSqr { };
99 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
100 V u[DESCRIPTOR::d];
101 MomentaF().computeU(cell, u);
102 //const V rho = MomentaF().computeRho(cell);
103 const V rho = cell.template getField<descriptors::RHO>();
104 const auto nablaRho = cell.template getField<descriptors::NABLARHO>();
105 V nablaRhoSqr = util::normSqr<V,DESCRIPTOR::d>(nablaRho);
106 V nablaRhoU = nablaRho[0]*u[0]+nablaRho[1]*u[1];
107 const V omega = parameters.template get<descriptors::OMEGA>();
108 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
109 V piNeqNorm = util::sqrt(piNeqNormSqr);
110 V preFactor = smagorinsky*smagorinsky
111 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
112 * 2 * util::sqrt(2);
114 V tauMol = V{1} / omega;
116 V tauTurb = V{0.5} * (util::sqrt(tauMol*tauMol + preFactor / rho * piNeqNorm) - tauMol);
118 V tauEff = tauMol + tauTurb;
119 return V{1} / tauEff;
120 }
121
122 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
123 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
124 parameters.template set<descriptors::OMEGA>(
125 computeEffectiveOmega(cell, parameters));
126 return CollisionO().apply(cell, parameters);
127 }
128};
129
130template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
131struct LocalSmagorinskyEffectiveOmega {
132 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
133 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
134
135 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
136 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
137 V piNeqNormSqr { };
138 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
139 const V rho = MomentaF().computeRho(cell);
140 const V omega = parameters.template get<descriptors::OMEGA>();
141 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
142 V piNeqNorm = util::sqrt(piNeqNormSqr);
143 V preFactor = smagorinsky*smagorinsky
144 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
145 * 2 * util::sqrt(2);
147 V tauMol = V{1} / omega;
149 V tauTurb = V{0.5} * (util::sqrt(tauMol*tauMol + preFactor / rho * piNeqNorm) - tauMol);
151 V tauEff = tauMol + tauTurb;
152 cell.template setField<descriptors::OMEGA>(V{1} / tauEff);
153 return V{1} / tauEff;
154 }
155
156 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
157 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
158 parameters.template set<descriptors::OMEGA>(
159 computeEffectiveOmega(cell, parameters));
160 return CollisionO().apply(cell, parameters);
161 }
162};
163
164template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
165struct LocalVanDriestSmagorinskyEffectiveOmega {
166 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
167 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
168
169 static constexpr bool is_vectorizable = false;
170
171 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
172 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
173 V piNeqNormSqr { };
174 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
175 const V rho = MomentaF().computeRho(cell);
176 const V omega = parameters.template get<descriptors::OMEGA>();
177 V tauEff {};
178 V nuT = cell.template getField<descriptors::VISCOSITY>();
179 if ( nuT != V(0) ) {
180 tauEff = nuT*descriptors::invCs2<V,DESCRIPTOR>() + V(0.5);
181 } else {
182 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
183 V piNeqNorm = util::sqrt(piNeqNormSqr);
184 V preFactor = smagorinsky*smagorinsky
185 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
186 * 2 * util::sqrt(2);
188 V tauMol = V{1} / omega;
190 V tauTurb = V{0.5} * (util::sqrt(tauMol*tauMol + preFactor / rho * piNeqNorm) - tauMol);
192 tauEff = tauMol + tauTurb;
193 }
194 cell.template setField<descriptors::OMEGA>(V{1} / tauEff);
195 return V{1} / tauEff;
196 }
197
198 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
199 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
200 parameters.template set<descriptors::OMEGA>(
201 computeEffectiveOmega(cell, parameters));
202 return CollisionO().apply(cell, parameters);
203 }
204};
205
206template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
207struct ShearSmagorinskyEffectiveOmega {
208 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
209 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
210
211 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
212 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
213 V piNeqNormSqr { };
214 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
215 const V rho = MomentaF().computeRho(cell);
216 const auto iT = parameters.template get<descriptors::LATTICE_TIME>();
217 const V omega = parameters.template get<descriptors::OMEGA>();
218 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
219 V avShear = cell.template getField<descriptors::AV_SHEAR>();
220 V piNeqNorm = util::sqrt(piNeqNormSqr);
221 V preFactor = smagorinsky*smagorinsky
222 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
223 * 2 * util::sqrt(2);
224 avShear = (avShear*iT + piNeqNorm) / (iT+1);
225 V tauMol = V{1} / omega;
226 V piNeqNormSISM = piNeqNorm - avShear;
227 V tauTurb = V{0.5} * (util::sqrt(tauMol*tauMol+(preFactor*piNeqNormSISM/rho))-tauMol);
228 V tauEff = tauMol + tauTurb;
229 return V{1} / tauEff;
230 }
231
232 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
233 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
234 parameters.template set<descriptors::OMEGA>(
235 computeEffectiveOmega(cell, parameters));
236 const auto iT = parameters.template get<descriptors::LATTICE_TIME>();
237 V piNeqNormSqr { };
238 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
239 V avShear = cell.template getField<descriptors::AV_SHEAR>();
240 avShear = (avShear*iT + util::sqrt(piNeqNormSqr)) / (iT+1);
241 cell.template setField<descriptors::AV_SHEAR>(avShear);
242 return CollisionO().apply(cell, parameters);
243 }
244};
245
246template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
247struct ConStrainSmagorinskyEffectiveOmega {
248 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
249 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
250
251 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
252 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
253 V pi[util::TensorVal<DESCRIPTOR>::n] { };
254 MomentaF().computeStress(cell, pi);
255 V piNeqNormSqr { };
256 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
257 const V rho = MomentaF().computeRho(cell);
258 const V omega = parameters.template get<descriptors::OMEGA>();
259 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
260 V piNeqNorm = util::sqrt(2*piNeqNormSqr);
261 V S[util::TensorVal<DESCRIPTOR>::n];
262 V cs2 = V{1} / descriptors::invCs2<V,DESCRIPTOR>();
263 V tauMol = V{1} / omega;
264 //Strain Tensor
265 V Phi = (-0.5*(-rho*tauMol*cs2+util::sqrt(rho*rho*tauMol*tauMol*cs2*cs2+V{2}*(smagorinsky*smagorinsky)*rho*piNeqNorm))/(smagorinsky*smagorinsky*rho*piNeqNorm));
266 for (int n=0; n < util::TensorVal<DESCRIPTOR>::n; ++n) {
267 S[n] = Phi*pi[n];
268 }
269 //Strain Tensor Norm
270 V SNormSqr = S[0]*S[0] + 2.0*S[1]*S[1] + S[2]*S[2];
271 if constexpr (util::TensorVal<DESCRIPTOR>::n == 6) {
272 SNormSqr += S[2]*S[2] + S[3]*S[3] + 2.0*S[4]*S[4] + S[5]*S[5];
273 }
274 V SNorm = util::sqrt(2*SNormSqr);
276 V tauTurb = smagorinsky*smagorinsky*SNorm/cs2;
278 V tauEff = tauMol + tauTurb;
279 return V{1} / tauEff;
280 }
281
282 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
283 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform{
284 parameters.template set<descriptors::OMEGA>(
285 computeEffectiveOmega(cell, parameters));
286 return CollisionO().apply(cell, parameters);
287 }
288
289};
290
291template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
292struct ConSmagorinskyEffectiveOmega {
293 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
294 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
295
296 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
297 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
298 V pi[util::TensorVal<DESCRIPTOR>::n] { };
299 MomentaF().computeStress(cell, pi);
300 V piNeqNormSqr { };
301 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
302 const V rho = MomentaF().computeRho(cell);
303 const V omega = parameters.template get<descriptors::OMEGA>();
304 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
305 V H[util::TensorVal<DESCRIPTOR >::n];
306 V conSmagoR[DESCRIPTOR::q];
307 V S[util::TensorVal<DESCRIPTOR >::n];
308 V tauMol = V{1} / omega;
309 V cs2 = V{1} / descriptors::invCs2<V,DESCRIPTOR>();
310 V piNeqNorm = util::sqrt(2*piNeqNormSqr);
311 V Phi = (-0.5*(-rho*tauMol*cs2+util::sqrt(rho*rho*tauMol*tauMol*cs2*cs2+V{2}*(smagorinsky*smagorinsky)*rho*piNeqNorm))/(smagorinsky*smagorinsky*rho*piNeqNorm));
312 for (int n=0; n < util::TensorVal<DESCRIPTOR>::n; ++n) {
313 S[n] = Phi*pi[n];
314 }
315 //Strain rate Tensor Norm
316 V SNormSqr = S[0]*S[0] + 2.0*S[1]*S[1] + S[2]*S[2];
317 if constexpr (util::TensorVal<DESCRIPTOR>::n == 6) {
318 SNormSqr += S[2]*S[2] + S[3]*S[3] + 2.0*S[4]*S[4] + S[5]*S[5];
319 }
320 V SNorm = util::sqrt(2*SNormSqr);
321 V preFactor = smagorinsky*smagorinsky
322 * descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
323 * 2 * util::sqrt(2);
324 //consistent Samagorinsky additional R term
325 //for (int q=0; q < DESCRIPTOR::q; ++q) {
326 {
327 unsigned q = 0;
328 V t = descriptors::t<V,DESCRIPTOR>(q);
329 //Hermite-Polynom H = c*c-cs^2*kronDelta
330 H[0] = descriptors::c<DESCRIPTOR>(q,0)*descriptors::c<DESCRIPTOR>(q,0)-cs2;
331 H[1] = descriptors::c<DESCRIPTOR>(q,0)*descriptors::c<DESCRIPTOR>(q,1);
332 H[2] = descriptors::c<DESCRIPTOR>(q,1)*descriptors::c<DESCRIPTOR>(q,1)-cs2;
333 if constexpr (util::TensorVal<DESCRIPTOR>::n == 6) {
334 H[2] = descriptors::c<DESCRIPTOR>(q,0)*descriptors::c<DESCRIPTOR>(q,2);
335 H[3] = descriptors::c<DESCRIPTOR>(q,1)*descriptors::c<DESCRIPTOR>(q,1)-cs2;
336 H[4] = descriptors::c<DESCRIPTOR>(q,1)*descriptors::c<DESCRIPTOR>(q,2);
337 H[5] = descriptors::c<DESCRIPTOR>(q,2)*descriptors::c<DESCRIPTOR>(q,2)-cs2;
338 }
339 //contraction or scalar product H*S
340 V contractHS = H[0]*S[0] + 2.0*H[1]*S[1] + H[2]*S[2];
341 if constexpr (util::TensorVal<DESCRIPTOR>::n == 6) {
342 contractHS += H[2]*S[2] + H[3]*S[3] + 2.0*H[4]*S[4] + H[5]*S[5];
343 }
344 //additional term
345 conSmagoR[q] = t*preFactor*SNorm*contractHS;
346 }
347 return conSmagoR[0];
348 }
349
350 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
351 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
352 parameters.template set<descriptors::OMEGA>(
353 computeEffectiveOmega(cell, parameters));
354 return CollisionO().apply(cell, parameters);
355 }
356};
357
358template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
359struct WaleEffectiveOmega {
360 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
361 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
362
363 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
364 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
365 V piNeqNormSqr { };
366 MomentaF().computePiNeqNormSqr(cell, piNeqNormSqr);
367 const V omega = parameters.template get<descriptors::OMEGA>();
368 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
369 const auto velocityGradient = cell.template getField<descriptors::VELO_GRAD>();
370 V preFactor = smagorinsky*smagorinsky;
371 // velocity gradient tensor
372 V g[3][3];
373 for (unsigned i=0; i < 3; i++) {
374 for (unsigned j=0; j < 3; j++) {
375 g[i][j] = velocityGradient[i*3 + j];
376 }
377 }
378 // strain rate tensor
379 V s[3][3];
380 for (unsigned i = 0; i < 3; i++) {
381 for (unsigned j = 0; j < 3; j++) {
382 s[i][j] = (g[i][j] + g[j][i]) / V{2};
383 }
384 }
385 // traceless symmetric part of the square of the velocity gradient tensor
386 V G[3][3];
387 for (unsigned i = 0; i < 3; i++) {
388 for (unsigned j = 0; j < 3; j++) {
389 G[i][j] = V{0};
390 }
391 }
392 for (unsigned i = 0; i < 3; i++) {
393 for (unsigned j = 0; j < 3; j++) {
394 for (unsigned k = 0; k < 3; k++) {
395 G[i][j] += (g[i][k]*g[k][j] + g[j][k]*g[k][i]) / V{2};
396 }
397 }
398 }
399 V trace{};
400 for (unsigned i = 0; i < 3; i++) {
401 trace += V{1}/V{3} * g[i][i] * g[i][i];
402 }
403 for (unsigned i = 0; i < 3; i++) {
404 G[i][i] -= trace;
405 }
406 // inner product of the traceless symmetric part of the square of the velocity gradient tensor
407 V G_ip{};
408 for (unsigned i = 0; i < 3; i++) {
409 for (unsigned j = 0; j < 3; j++) {
410 G_ip = G[i][j] * G[i][j];
411 }
412 }
413 // inner product of the strain rate
414 V s_ip{};
415 for (unsigned i = 0; i < 3; i++) {
416 for (unsigned j = 0; j < 3; j++) {
417 s_ip = s[i][j] * s[i][j];
418 }
419 }
420 // Turbulent relaxation time
421 V tauTurb = V{3} * preFactor * (util::pow(G_ip,1.5) / (util::pow(s_ip,2.5) + util::pow(G_ip,1.25)));
422 if ((util::pow(s_ip,2.5) + util::pow(G_ip,1.25)) == 0) {
423 tauTurb = 0;
424 }
425 // Physical turbulent viscosity must be equal or higher that zero
426 if (tauTurb < 0) {
427 tauTurb = 0;
428 }
430 V tauMol = V{1} / omega;
432 V tauEff = tauMol + tauTurb;
433 return V{1} / tauEff;
434 }
435
436 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
437 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
438 parameters.template set<descriptors::OMEGA>(
439 computeEffectiveOmega(cell, parameters));
440 return CollisionO().apply(cell, parameters);
441 }
442};
443
444template <typename COLLISION, typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
445struct KrauseEffectiveOmega {
446 using MomentaF = typename MOMENTA::template type<DESCRIPTOR>;
447 using CollisionO = typename COLLISION::template type<DESCRIPTOR, MOMENTA, EQUILIBRIUM>;
448
449 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
450 V computeEffectiveOmega(CELL& cell, PARAMETERS& parameters) any_platform {
451 V rho, u[DESCRIPTOR::d], fNeq[DESCRIPTOR::q] { };
452 MomentaF().computeRhoU(cell, rho, u);
453 lbm<DESCRIPTOR>::computeFneq(cell, fNeq, rho, u);
454 const V omega = parameters.template get<descriptors::OMEGA>();
455 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
456 const V preFactor = smagorinsky*smagorinsky
457 * 3*descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
458 * 2*util::sqrt(2);
459 FieldD<V,DESCRIPTOR,typename COLLISION::OMEGA> omegaEff;
460 const V tauMol = V{1} / omega;
461 for (unsigned iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
462 V tauTurb = V{0.5}*(util::sqrt(tauMol*tauMol + preFactor/rho * util::fabs(fNeq[iPop])) - tauMol);
463 omegaEff[iPop] = V{1} / (tauMol + tauTurb);
464 }
465 V avgOmegaEff = 0.;
466 for (int iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
467 avgOmegaEff += omegaEff[iPop];
468 }
469 avgOmegaEff /= DESCRIPTOR::q;
470 return avgOmegaEff;
471 }
472
473 template <concepts::Cell CELL, concepts::Parameters PARAMETERS, typename V=typename CELL::value_t>
474 CellStatistic<V> apply(CELL& cell, PARAMETERS& parameters) any_platform {
475 V rho, u[DESCRIPTOR::d], fNeq[DESCRIPTOR::q] { };
476 MomentaF().computeRhoU(cell, rho, u);
477 lbm<DESCRIPTOR>::computeFneq(cell, fNeq, rho, u); // TODO: Use EQUILIBRIUM instead of 2nd order
478 const V omega = parameters.template get<descriptors::OMEGA>();
479 const V smagorinsky = parameters.template get<collision::LES::SMAGORINSKY>();
480 const V preFactor = smagorinsky*smagorinsky
481 * 3*descriptors::invCs2<V,DESCRIPTOR>()*descriptors::invCs2<V,DESCRIPTOR>()
482 * 2*util::sqrt(2);
483 FieldD<V,DESCRIPTOR,typename COLLISION::OMEGA> omegaEff;
484 const V tauMol = V{1} / omega;
485 for (unsigned iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
486 V tauTurb = V{0.5}*(util::sqrt(tauMol*tauMol + preFactor/rho * util::fabs(fNeq[iPop])) - tauMol);
487 omegaEff[iPop] = V{1} / (tauMol + tauTurb);
488 }
489 parameters.template set<typename COLLISION::OMEGA>(omegaEff);
490 return CollisionO().apply(cell, parameters);
491 }
492};
493
494
495}
496
498template <typename COLLISION>
499struct SmagorinskyEffectiveOmega {
500 using parameters = typename COLLISION::parameters::template include<
501 descriptors::OMEGA, LES::SMAGORINSKY
502 >;
503
504 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
505 "COLLISION must be parametrized using relaxation frequency OMEGA");
506
507 static std::string getName() {
508 return "SmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
509 }
510
511 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
512 using type = detail::SmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
513};
514
516template <typename COLLISION>
517struct IncompressibleSmagorinskyEffectiveOmega {
518 using parameters = typename COLLISION::parameters::template include<
519 descriptors::OMEGA, LES::SMAGORINSKY
520 >;
521
522 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
523 "COLLISION must be parametrized using relaxation frequency OMEGA");
524
525 static std::string getName() {
526 return "SmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
527 }
528
529 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
530 using type = detail::IncompressibleSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
531};
532
533template <typename COLLISION>
534struct LocalSmagorinskyEffectiveOmega {
535 using parameters = typename COLLISION::parameters::template include<
536 descriptors::OMEGA, LES::SMAGORINSKY
537 >;
538
539 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
540 "COLLISION must be parametrized using relaxation frequency OMEGA");
541
542 static std::string getName() {
543 return "LocalSmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
544 }
545
546 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
547 using type = detail::LocalSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
548};
549
550template <typename COLLISION>
551struct LocalVanDriestSmagorinskyEffectiveOmega {
552 using parameters = typename COLLISION::parameters::template include<
553 descriptors::OMEGA, LES::SMAGORINSKY
554 >;
555
556 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
557 "COLLISION must be parametrized using relaxation frequency OMEGA");
558
559 static std::string getName() {
560 return "LocalVanDriestSmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
561 }
562
563 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
564 using type = detail::LocalVanDriestSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
565};
566
568template <typename COLLISION>
569struct ShearSmagorinskyEffectiveOmega {
570 using parameters = typename COLLISION::parameters::template include<
571 descriptors::LATTICE_TIME, descriptors::OMEGA, LES::SMAGORINSKY
572 >;
573
574 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
575 "COLLISION must be parametrized using relaxation frequency OMEGA");
576
577 static std::string getName() {
578 return "ShearSmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
579 }
580
581 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
582 using type = detail::ShearSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
583};
584
586template <typename COLLISION>
587struct ConStrainSmagorinskyEffectiveOmega {
588 using parameters = typename COLLISION::parameters::template include<
589 descriptors::OMEGA, LES::SMAGORINSKY
590 >;
591
592 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
593 "COLLISION must be parametrized using relaxation frequency OMEGA");
594
595 static std::string getName() {
596 return "ConStrainSmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
597 }
598
599 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
600 using type = detail::ConStrainSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
601};
602
604template <typename COLLISION>
605struct ConSmagorinskyEffectiveOmega {
606 using parameters = typename COLLISION::parameters::template include<
607 descriptors::OMEGA, LES::SMAGORINSKY
608 >;
609
610 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
611 "COLLISION must be parametrized using relaxation frequency OMEGA");
612
613 static std::string getName() {
614 return "ConSmagorinskyEffectiveOmega<" + COLLISION::getName() + ">";
615 }
616
617 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
618 using type = detail::ConSmagorinskyEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
619};
620
622template <typename COLLISION>
623struct WaleEffectiveOmega {
624 using parameters = typename COLLISION::parameters::template include<
625 descriptors::OMEGA, LES::SMAGORINSKY
626 >;
627
628 static_assert(COLLISION::parameters::template contains<descriptors::OMEGA>(),
629 "COLLISION must be parametrized using relaxation frequency OMEGA");
630
631 static std::string getName() {
632 return "WaleEffectiveOmega<" + COLLISION::getName() + ">";
633 }
634
635 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
636 using type = detail::WaleEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
637};
638
640template <typename COLLISION>
641struct KrauseEffectiveOmega {
642 using parameters = typename COLLISION::parameters::template include<
643 descriptors::OMEGA, LES::SMAGORINSKY
644 >;
645
646 static std::string getName() {
647 return "KrauseEffectiveOmega<" + COLLISION::getName() + ">";
648 }
649
650 template <typename DESCRIPTOR, typename MOMENTA, typename EQUILIBRIUM>
651 using type = detail::KrauseEffectiveOmega<COLLISION,DESCRIPTOR,MOMENTA,EQUILIBRIUM>;
652};
653
654}
655
656}
657
658#endif
olb::Vector
Plain old scalar vector.
Definition vectorHelpers.h:41
olb::descriptors::invCs2
constexpr T invCs2() any_platform
Definition functions.h:107
olb::descriptors::t
constexpr T t(unsigned iPop, tag::CUM) any_platform
Definition cum.h:108
olb::descriptors::c
constexpr int c(unsigned iPop, unsigned iDim) any_platform
Definition functions.h:83
olb::util::sqrt
Expr sqrt(Expr x)
Definition expr.cpp:225
olb::util::normSqr
auto normSqr(const ARRAY_LIKE &u) any_platform
Compute norm square of a d-dimensional vector.
Definition util.h:145
olb::util::pow
Expr pow(Expr base, Expr exp)
Definition expr.cpp:235
olb::util::fabs
Expr fabs(Expr x)
Definition expr.cpp:230
olb
Top level namespace for all of OpenLB.
Definition advectionDiffusionDirichlet.h:27
olb::FieldD
Vector< typename FIELD::template value_type< T >, DESCRIPTOR::template size< FIELD >() > FieldD
Vector storing a single field instance.
Definition vector.h:480
any_platform
#define any_platform
Define preprocessor macros for device-side functions, constant storage.
Definition platform.h:77
olb::CellStatistic
Return value of any collision.
Definition interface.h:45
olb::collision::ConSmagorinskyEffectiveOmega
Compute dynamics parameter OMEGA locally using Consistent Smagorinsky LES model.
Definition collisionLES.h:605
olb::collision::ConSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:606
olb::collision::ConStrainSmagorinskyEffectiveOmega
Compute dynamics parameter OMEGA locally using Consistent Strain Smagorinsky LES model.
Definition collisionLES.h:587
olb::collision::ConStrainSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:588
olb::collision::IncompressibleSmagorinskyEffectiveOmega
Compute dynamics parameter OMEGA locally using Smagorinsky LES model.
Definition collisionLES.h:517
olb::collision::IncompressibleSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:518
olb::collision::KrauseEffectiveOmega
Compute dynamics parameter OMEGA locally using Krause LES model.
Definition collisionLES.h:641
olb::collision::KrauseEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:642
olb::collision::LES::SMAGORINSKY::isValid
static constexpr auto isValid(FieldD< T, DESCRIPTOR, FIELD > value)
Definition collisionLES.h:44
olb::collision::LES::SMAGORINSKY::getInitialValue
static constexpr auto getInitialValue()
Definition collisionLES.h:39
olb::collision::LocalSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:535
olb::collision::LocalVanDriestSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:552
olb::collision::ShearSmagorinskyEffectiveOmega
Compute dynamics parameter OMEGA locally using Shear Smagorinsky LES model.
Definition collisionLES.h:569
olb::collision::ShearSmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::LATTICE_TIME, descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:570
olb::collision::SmagorinskyEffectiveOmega
Compute dynamics parameter OMEGA locally using Smagorinsky LES model.
Definition collisionLES.h:499
olb::collision::SmagorinskyEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:500
olb::collision::WaleEffectiveOmega
Compute dynamics parameter OMEGA locally using WALE.
Definition collisionLES.h:623
olb::collision::WaleEffectiveOmega::parameters
typename COLLISION::parameters::template include< descriptors::OMEGA, LES::SMAGORINSKY > parameters
Definition collisionLES.h:624
olb::collision::detail::ConSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:294
olb::collision::detail::ConSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:293
olb::collision::detail::ConSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:297
olb::collision::detail::ConSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:351
olb::collision::detail::ConStrainSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:252
olb::collision::detail::ConStrainSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:248
olb::collision::detail::ConStrainSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:249
olb::collision::detail::ConStrainSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:283
olb::collision::detail::IncompressibleSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:93
olb::collision::detail::IncompressibleSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:94
olb::collision::detail::IncompressibleSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:123
olb::collision::detail::IncompressibleSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:97
olb::collision::detail::KrauseEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:450
olb::collision::detail::KrauseEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:447
olb::collision::detail::KrauseEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:446
olb::collision::detail::KrauseEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:474
olb::collision::detail::LocalSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:157
olb::collision::detail::LocalSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:136
olb::collision::detail::LocalSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:132
olb::collision::detail::LocalSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:133
olb::collision::detail::LocalVanDriestSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:166
olb::collision::detail::LocalVanDriestSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:172
olb::collision::detail::LocalVanDriestSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:167
olb::collision::detail::LocalVanDriestSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:199
olb::collision::detail::ShearSmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:208
olb::collision::detail::ShearSmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:212
olb::collision::detail::ShearSmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:233
olb::collision::detail::ShearSmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:209
olb::collision::detail::SmagorinskyEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:64
olb::collision::detail::SmagorinskyEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:84
olb::collision::detail::SmagorinskyEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:60
olb::collision::detail::SmagorinskyEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:61
olb::collision::detail::WaleEffectiveOmega::apply
CellStatistic< V > apply(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:437
olb::collision::detail::WaleEffectiveOmega::CollisionO
typename COLLISION::template type< DESCRIPTOR, MOMENTA, EQUILIBRIUM > CollisionO
Definition collisionLES.h:361
olb::collision::detail::WaleEffectiveOmega::MomentaF
typename MOMENTA::template type< DESCRIPTOR > MomentaF
Definition collisionLES.h:360
olb::collision::detail::WaleEffectiveOmega::computeEffectiveOmega
V computeEffectiveOmega(CELL &cell, PARAMETERS &parameters) any_platform
Definition collisionLES.h:364
olb::descriptors::FIELD_BASE
Base of a field whose size is defined by [C_0,C_1,C_2]^T * [1,D,Q].
Definition fields.h:52
olb::descriptors::LATTICE_TIME
Definition fields.h:308
olb::lbm::computeFneq
static void computeFneq(CELL &cell, FNEQ &fNeq, const RHO &rho, const U &u) any_platform
Computation of non-equilibrium distribution.
Definition lbm.h:353
olb::util::TensorVal
Compute number of elements of a symmetric d-dimensional tensor.
Definition util.h:216
Generated on Sat May 17 2025 10:48:00 for OpenLB by doxygen 1.12.0