24#ifndef TURBULENT_F_3D_HH
25#define TURBULENT_F_3D_HH
33template <
typename T,
typename DESCRIPTOR>
38 _superGeometry(superGeometry), _indicator(indicator), _material(material)
43template <
typename T,
typename DESCRIPTOR>
46 int globIC = input[0];
53 std::vector<T> normalTemp(3,T());
54 std::vector<T> normal(3,T());
57 if (_superGeometry.get(input) == 1) {
58 for (
int iPop = 1; iPop < DESCRIPTOR::q; ++iPop) {
59 if (_superGeometry.get(input[0],
60 input[1] + descriptors::c<DESCRIPTOR>(iPop,0),
61 input[2] + descriptors::c<DESCRIPTOR>(iPop,1),
62 input[3] + descriptors::c<DESCRIPTOR>(iPop,2)) == _material) {
64 normalTemp[0] += descriptors::c<DESCRIPTOR>(iPop,0);
65 normalTemp[1] += descriptors::c<DESCRIPTOR>(iPop,1);
66 normalTemp[2] += descriptors::c<DESCRIPTOR>(iPop,2);
72 std::vector<T> physR(3, T());
73 _superGeometry.getCuboidGeometry().getPhysR(&(physR[0]), &(input[0]));
75 T voxelSize = _superGeometry.getCuboidGeometry().get(globIC).getDeltaR();
79 std::vector<T> direction(normal);
80 direction[0] = voxelSize*normal[0]*2.;
81 direction[1] = voxelSize*normal[1]*2.;
82 direction[2] = voxelSize*normal[2]*2.;
85 if ( _indicator.distance(distance, physR, direction) ) {
90 auto cell = this->_sLattice.getBlock(this->_sLattice.getLoadBalancer().loc(globIC)).get(locix, lociy, lociz);
91 cell.computeAllMomenta(rho, u, pi);
94 T omega = 1. / this->_converter.getLatticeRelaxationTime();
95 T dt = this->_converter.getConversionFactorTime();
96 T physFactor = -omega*descriptors::invCs2<T,DESCRIPTOR>()/rho/2./dt*this->_converter.getPhysDensity(rho)*this->_converter.getPhysViscosity();
99 T Rx = pi[0]*physFactor*normal[0] + pi[1]*physFactor*normal[1] + pi[2]*physFactor*normal[2];
100 T Ry = pi[1]*physFactor*normal[0] + pi[3]*physFactor*normal[1] + pi[4]*physFactor*normal[2];
101 T Rz = pi[2]*physFactor*normal[0] + pi[4]*physFactor*normal[1] + pi[5]*physFactor*normal[2];
104 T R_res_pressure = normal[0]*pi[0]*physFactor*normal[0] + normal[0]*pi[1]*physFactor*normal[1] + normal[0]*pi[2]*physFactor*normal[2]
105 +normal[1]*pi[1]*physFactor*normal[0] + normal[1]*pi[3]*physFactor*normal[1] + normal[1]*pi[4]*physFactor*normal[2]
106 +normal[2]*pi[2]*physFactor*normal[0] + normal[2]*pi[4]*physFactor*normal[1] + normal[2]*pi[5]*physFactor*normal[2];
108 Rx -= R_res_pressure *normal[0];
109 Ry -= R_res_pressure *normal[1];
110 Rz -= R_res_pressure *normal[2];
113 T u_tau =
util::sqrt(tau_wall/this->_converter.getPhysDensity(rho));
115 output[0] = distance*u_tau / this->_converter.getPhysViscosity();
124template <
typename T,
typename DESCRIPTOR>
127 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockFinDiff(blockFinDiff)
129 this->
getName() =
"VelocityGradientFD";
132template <
typename T,
typename DESCRIPTOR>
138 _blockFinDiff(output,input);
143template <
typename T,
typename DESCRIPTOR>
146 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockFinDiff(blockFinDiff)
148 this->
getName() =
"externalVelocityGradientFD";
151template <
typename T,
typename DESCRIPTOR>
157 _blockFinDiff(output,input);
162template <
typename T,
typename DESCRIPTOR>
165 _sVelocity(sLattice), _sFinDiff(sGeometry, _sVelocity, matNumber)
167 this->
getName() =
"VelocityGradientFD";
170 for (
int iC = 0; iC < maxC; iC++) {
176template <
typename T,
typename DESCRIPTOR>
179 _sVelocity(sLattice), _sFinDiff(sGeometry, _sVelocity, matNumber)
181 this->
getName() =
"externalVelocityGradientFD";
184 for (
int iC = 0; iC < maxC; iC++) {
190template <
typename T,
typename DESCRIPTOR>
193 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockFinDiff(blockFinDiff), _converter(converter)
195 this->
getName() =
"PhysVelocityGradientFD";
198template <
typename T,
typename DESCRIPTOR>
201 _blockFinDiff(output,input);
206template <
typename T,
typename DESCRIPTOR>
209 _sVelocity(sLattice, converter), _sFinDiff(sGeometry, _sVelocity, matNumber, converter), _converter(converter)
211 this->
getName() =
"PhysVelocityGradientFD";
214 for (
int iC = 0; iC < maxC; iC++) {
220template <
typename T,
typename DESCRIPTOR>
223 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockVeloGrad(blockVeloGrad)
225 this->
getName() =
"StrainRateFD";
228template <
typename T,
typename DESCRIPTOR>
232 _blockVeloGrad(velograd,input);
233 output[0] = velograd[0];
234 output[1] = 0.5 * (velograd[1] + velograd[3]);
235 output[2] = 0.5 * (velograd[2] + velograd[6]);
236 output[3] = output[1];
237 output[4] = velograd[4];
238 output[5] = 0.5 * (velograd[5] + velograd[7]);
239 output[6] = output[2];
240 output[7] = output[5];
241 output[8] = velograd[8];
246template <
typename T,
typename DESCRIPTOR>
249 _sVeloGrad(sGeometry, sLattice, matNumber)
251 this->
getName() =
"StrainRateFD";
254 for (
int iC = 0; iC < maxC; iC++) {
260template <
typename T,
typename DESCRIPTOR>
263 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockVeloGrad(blockVeloGrad), _converter(converter)
265 this->
getName() =
"PhysStrainRateFD";
268template <
typename T,
typename DESCRIPTOR>
272 _blockVeloGrad(velograd,input);
273 output[0] = velograd[0];
274 output[1] = 0.5 * (velograd[1] + velograd[3]);
275 output[2] = 0.5 * (velograd[2] + velograd[6]);
276 output[3] = output[1];
277 output[4] = velograd[4];
278 output[5] = 0.5 * (velograd[5] + velograd[7]);
279 output[6] = output[2];
280 output[7] = output[5];
281 output[8] = velograd[8];
287template <
typename T,
typename DESCRIPTOR>
290 _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
292 this->
getName() =
"PhysStrainRateFD";
295 for (
int iC = 0; iC < maxC; iC++) {
301template <
typename T,
typename DESCRIPTOR>
304 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 1), _blockVeloGrad(blockVeloGrad), _converter(converter)
306 this->
getName() =
"DissipationFD";
309template <
typename T,
typename DESCRIPTOR>
313 _blockVeloGrad(velograd,input);
314 output[0] = velograd[0] * velograd[0] + velograd[1] * velograd[1] + velograd[2] * velograd[2] +
315 velograd[3] * velograd[3] + velograd[4] * velograd[4] + velograd[5] * velograd[5] +
316 velograd[6] * velograd[6] + velograd[7] * velograd[7] + velograd[8] * velograd[8];
317 output[0] *= _converter.getLatticeViscosity();
323template <
typename T,
typename DESCRIPTOR>
326 _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
328 this->
getName() =
"DissipationFD";
331 for (
int iC = 0; iC < maxC; iC++) {
337template <
typename T,
typename DESCRIPTOR>
340 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 1), _blockVeloGrad(blockVeloGrad), _converter(converter)
342 this->
getName() =
"PhysDissipationFD";
345template <
typename T,
typename DESCRIPTOR>
349 _blockVeloGrad(velograd,input);
350 output[0] = velograd[0] * velograd[0] + velograd[1] * velograd[1] + velograd[2] * velograd[2] +
351 velograd[3] * velograd[3] + velograd[4] * velograd[4] + velograd[5] * velograd[5] +
352 velograd[6] * velograd[6] + velograd[7] * velograd[7] + velograd[8] * velograd[8];
353 output[0] *= _converter.getPhysViscosity();
359template <
typename T,
typename DESCRIPTOR>
362 _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
364 this->
getName() =
"PhysDissipationFD";
367 for (
int iC = 0; iC < maxC; iC++) {
373template <
typename T,
typename DESCRIPTOR>
380 _blockVeloGrad(blockVeloGrad),
381 _converter(converter),
382 _effectiveOmegaF(effectiveOmegaF)
384 this->
getName() =
"PhysEffectiveDissipationFD";
387template <
typename T,
typename DESCRIPTOR>
391 _blockVeloGrad(velograd,input);
392 output[0] = velograd[0] * velograd[0] + velograd[1] * velograd[1] + velograd[2] * velograd[2] +
393 velograd[3] * velograd[3] + velograd[4] * velograd[4] + velograd[5] * velograd[5] +
394 velograd[6] * velograd[6] + velograd[7] * velograd[7] + velograd[8] * velograd[8];
396 auto cell = this->_blockLattice.get(input[0], input[1], input[2]);
397 T omegaEff = _effectiveOmegaF(cell);
398 T nuEff = ((1./omegaEff)-0.5)/descriptors::invCs2<T,DESCRIPTOR>();
399 output[0] *= _converter.getPhysViscosity( nuEff );
405template <
typename T,
typename DESCRIPTOR>
409 :
SuperLatticeF3D<T,DESCRIPTOR>(sLattice,1), _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
411 this->
getName() =
"PhysEffectiveDissipationFD";
414 for (
int iC = 0; iC < maxC; iC++) {
420template <
typename T,
typename DESCRIPTOR>
423 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 3), _blockVeloGrad(blockVeloGrad)
425 this->
getName() =
"VorticityFD";
428template <
typename T,
typename DESCRIPTOR>
432 _blockVeloGrad(velograd,input);
433 output[0] = velograd[7] - velograd[5];
434 output[1] = velograd[2] - velograd[6];
435 output[2] = velograd[3] - velograd[1];
440template <
typename T,
typename DESCRIPTOR>
443 _sVeloGrad(sGeometry, sLattice, matNumber)
445 this->
getName() =
"VorticityFD";
448 for (
int iC = 0; iC < maxC; iC++) {
454template <
typename T,
typename DESCRIPTOR>
457 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 3), _blockVeloGrad(blockVeloGrad), _converter(converter)
459 this->
getName() =
"PhysVorticityFD";
462template <
typename T,
typename DESCRIPTOR>
466 _blockVeloGrad(velograd,input);
467 output[0] = velograd[7] - velograd[5];
468 output[1] = velograd[2] - velograd[6];
469 output[2] = velograd[3] - velograd[1];
474template <
typename T,
typename DESCRIPTOR>
477 _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
479 this->
getName() =
"PhysVorticityFD";
482 for (
int iC = 0; iC < maxC; iC++) {
488template <
typename T,
typename DESCRIPTOR>
491 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 9), _blockStrainRate(blockStrainRate), _converter(converter)
493 this->
getName() =
"PhysStressFD";
496template <
typename T,
typename DESCRIPTOR>
499 _blockStrainRate(output,input);
500 for (
int i = 0; i < 9; i++) {
501 output[i] /= _converter.getPhysViscosity() * _converter.getPhysDensity();
507template <
typename T,
typename DESCRIPTOR>
510 _sStrainRate(sGeometry, sLattice, matNumber, converter), _converter(converter)
512 this->
getName() =
"PhysStressFD";
515 for (
int iC = 0; iC < maxC; iC++) {
521template<
typename T,
typename DESCRIPTOR>
523 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 1), _blockVelocity(blockVelocity)
525 this->
getName() =
"IsotropicHomogeneousTKE";
528template<
typename T,
typename DESCRIPTOR>
532 T data[_blockVelocity.getTargetDim()];
533 _blockVelocity(data,input);
534 for (
int i = 0; i < _blockVelocity.getTargetDim(); ++i) {
535 output[0] += data[i] * data[i];
537 output[0] = 0.5 * output[0];
542template<
typename T,
typename DESCRIPTOR>
545 :
SuperLatticeF3D<T,DESCRIPTOR>(sLattice, 1), _sVelocity(sLattice, converter), _converter(converter)
548 this->
getName() =
"IsotropicHomogeneousTKE";
549 int maxC = this->
_sLattice.getLoadBalancer().size();
551 for (
int iC = 0; iC < maxC; iC++) {
557template <
typename T,
typename DESCRIPTOR>
560 :
BlockLatticeF3D<T,DESCRIPTOR>(blockLattice, 1), _blockVeloGrad(blockVeloGrad), _converter(converter)
562 this->
getName() =
"PhysEnstrophyFD";
565template <
typename T,
typename DESCRIPTOR>
569 _blockVeloGrad(velograd,input);
570 output[0] = 0.5 * (
util::pow(velograd[7] - velograd[5], 2) +
util::pow(velograd[2] - velograd[6], 2) +
util::pow(velograd[3] - velograd[1], 2) );
575template <
typename T,
typename DESCRIPTOR>
578 _sVeloGrad(sGeometry, sLattice, matNumber, converter), _converter(converter)
580 this->
getName() =
"PhysEnstrophyFD";
583 for (
int iC = 0; iC < maxC; iC++) {
represents all functors that operate on a cuboid in general, mother class of BlockLatticeF,...
functor that returns pointwise the turbulent, kinetic energy
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
BlockIsotropicHomogeneousTKE3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &f)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
BlockLatticeDissipationFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter)
BlockLatticeExternalVelocityGradientFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
represents all functors that operate on a DESCRIPTOR in general, e.g. getVelocity(),...
BlockLatticePhysDissipationFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &_converter)
bool operator()(T output[], const int input[]) override
has to be implemented for 'every' derived class
BlockLatticePhysEffectiveDissipationFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter, std::function< T(Cell< T, DESCRIPTOR > &)> effectiveOmegaF)
bool operator()(T output[], const int input[]) override
has to be implemented for 'every' derived class
functor that returns pointwise the enstrophy
BlockLatticePhysEnstrophyFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockVeloGrad, const UnitConverter< T, DESCRIPTOR > &converter)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
BlockLatticePhysStrainRateFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter)
BlockLatticePhysStressFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
bool operator()(T output[], const int input[]) override
has to be implemented for 'every' derived class
BlockLatticePhysVelocityGradientFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter)
BlockLatticePhysVorticityFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor, const UnitConverter< T, DESCRIPTOR > &converter)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
BlockLatticeStrainRateFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor)
functor to get pointwise explicit filtering on local lattice, if globIC is not on the local processor...
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
BlockLatticeVelocityGradientFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor)
BlockLatticeVorticityFD3D(BlockLattice< T, DESCRIPTOR > &blockLattice, BlockF3D< T > &blockFunctor)
bool operator()(T output[], const int input[])
has to be implemented for 'every' derived class
Platform-abstracted block lattice for external access and inter-block interaction.
Highest-level interface to Cell data.
std::string & getName()
read and write access to name
IndicatorF3D is an application from .
SuperStructure< T, 3 > & _superStructure
std::vector< std::unique_ptr< BlockF3D< T > > > _blockF
Super functors may consist of several BlockF3D<W> derived functors.
Representation of a statistic for a parallel 2D geometry.
SuperIsotropicHomogeneousTKE3D(SuperLattice< T, DESCRIPTOR > &sLattice, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticeDissipationFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticeExternalVelocityGradientFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber)
represents all functors that operate on a SuperLattice in general, e.g. getVelocity(),...
SuperLattice< T, DESCRIPTOR > & _sLattice
SuperLatticePhysDissipationFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticePhysEffectiveDissipationFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter, std::function< T(Cell< T, DESCRIPTOR > &)> effectiveOmegaF)
SuperLatticePhysEnstrophyFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
represents all functors that operate on a DESCRIPTOR with output in Phys, e.g. physVelocity(),...
SuperLatticePhysStrainRateFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticePhysStressFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticePhysVelocityGradientFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticePhysVorticityFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber, const UnitConverter< T, DESCRIPTOR > &converter)
SuperLatticeStrainRateFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber)
SuperLatticeVelocityGradientFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber)
SuperLatticeVorticityFD3D(SuperGeometry< T, 3 > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, std::list< int > &matNumber)
bool operator()(T output[], const int input[]) override
SuperLatticeYplus3D(SuperLattice< T, DESCRIPTOR > &sLattice, const UnitConverter< T, DESCRIPTOR > &converter, SuperGeometry< T, 3 > &superGeometry, IndicatorF3D< T > &indicator, const int material)
Super class maintaining block lattices for a cuboid decomposition.
LoadBalancer< T > & getLoadBalancer()
Read and write access to the load balancer.
Conversion between physical and lattice units, as well as discretization.
int getRank() const
Returns the process ID.
cpu::simd::Pack< T > sqrt(cpu::simd::Pack< T > value)
cpu::simd::Pack< T > pow(cpu::simd::Pack< T > base, cpu::simd::Pack< T > exp)
Vector< T, D > normalize(const Vector< T, D > &a)
Top level namespace for all of OpenLB.