OpenLB 1.8.1
Loading...
Searching...
No Matches
olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE > Class Template Referencefinal

#include <latticeSpheroidLiftForce.h>

+ Inheritance diagram for olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >:
+ Collaboration diagram for olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >:

Public Member Functions

 BlockLatticeSpheroidLiftForce (BlockLattice< T, DESCRIPTOR > &blockLattice, const BlockGeometry< T, DESCRIPTOR::d > &blockGeometry, particles::ParticleSystem< T, PARTICLETYPE > &particleSystem, const UnitConverter< T, DESCRIPTOR > &converter, PhysR< T, DESCRIPTOR::d > cellMin=PhysR< T, DESCRIPTOR::d >(0.), PhysR< T, DESCRIPTOR::d > cellMax=PhysR< T, DESCRIPTOR::d >(0.), Vector< bool, DESCRIPTOR::d > periodic=Vector< bool, DESCRIPTOR::d >(false), std::size_t iP0=0)
 
void evaluate (T output[], particles::Particle< T, PARTICLETYPE > &particle, int iP)
 
bool operator() (T output[], const int input[]) override
 

Static Public Attributes

static constexpr bool serializeForce = false
 

Detailed Description

template<typename T, typename DESCRIPTOR, typename PARTICLETYPE>
class olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >

Definition at line 62 of file latticeSpheroidLiftForce.h.

Constructor & Destructor Documentation

◆ BlockLatticeSpheroidLiftForce()

template<typename T , typename DESCRIPTOR , typename PARTICLETYPE >
olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >::BlockLatticeSpheroidLiftForce ( BlockLattice< T, DESCRIPTOR > & blockLattice,
const BlockGeometry< T, DESCRIPTOR::d > & blockGeometry,
particles::ParticleSystem< T, PARTICLETYPE > & particleSystem,
const UnitConverter< T, DESCRIPTOR > & converter,
PhysR< T, DESCRIPTOR::d > cellMin = PhysR<T,DESCRIPTOR::d> (0.),
PhysR< T, DESCRIPTOR::d > cellMax = PhysR<T,DESCRIPTOR::d> (0.),
Vector< bool, DESCRIPTOR::d > periodic = Vector<bool,DESCRIPTOR::d> (false),
std::size_t iP0 = 0 )

Definition at line 34 of file latticeSpheroidLiftForce.hh.

42 : BlockLatticePhysF<T,DESCRIPTOR>(blockLattice, converter,
43 (DESCRIPTOR::d)*(particleSystem.size()-iP0)),
44 _blockGeometry(blockGeometry), _blockLattice(blockLattice),
45 _particleSystem(particleSystem),
46 _cellMin(cellMin), _cellMax(cellMax), _periodic(periodic), _iP0(iP0)
47{
48 this->getName() = "physStokesSpheroidLiftForce";
49}
std::conditional_t< DESCRIPTOR::d==2, BlockLatticePhysF2D< T, DESCRIPTOR >, BlockLatticePhysF3D< T, DESCRIPTOR > > BlockLatticePhysF
Definition aliases.h:328
std::string getName(OperatorScope scope)
Returns human-readable name of scope.

References olb::getName().

+ Here is the call graph for this function:

Member Function Documentation

◆ evaluate()

template<typename T , typename DESCRIPTOR , typename PARTICLETYPE >
void olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >::evaluate ( T output[],
particles::Particle< T, PARTICLETYPE > & particle,
int iP )

Definition at line 52 of file latticeSpheroidLiftForce.hh.

53{
54 constexpr unsigned D = DESCRIPTOR::d;
55
56 using namespace descriptors;
57 using namespace eler;
58 static_assert(PARTICLETYPE::template providesNested<MOBILITY, FLUIDVEL>(), "Field MOBILITY:FLUIDVEL has to be provided");
59 //Retrieve particle quantities
60 Vector<T,D> position = particle.template getField<GENERAL,POSITION>();
61 Vector<T,D> velocity = particle.template getField<MOBILITY,VELOCITY>();
62 T* positionArray = position.data();
63
64 //TODO: check, whether creation can be avoided by using a second _blockF list
65 const auto& cuboid = _blockGeometry.getCuboid();
67 _blockLattice, this->_converter, cuboid);
68
69 //Calculate general constants
70
71 //Calculate particle coefficiants
72
73
74 T fluidVelArray[D] = {0.};
75
76 //Check whether inside cuboid (when not parallelized)
77 bool inside = true;
79 inside = cuboid.checkPoint(position);
80 }
81 if (inside){
82 //Force computation
83 //Retrieve interpolated velocity at position
84 blockInterpPhysVelF(fluidVelArray, positionArray);
85 particle.template setField<MOBILITY, FLUIDVEL>(fluidVelArray);
86 Vector<T,3> fvel = particle.template getField<MOBILITY, FLUIDVEL>();
87
88 Vector<T,9> rotmat = particle.template getField<NUMERICPROPERTIES,ROT_MATRIX>();
89
90 Vector<T,9> Khathat (0.);
91Khathat = eler::matrixMultiply(eler::inverseTransformMatrix(rotmat),
92eler::matrixMultiply(particle.template getField<NUMERICPROPERTIES,TRANSLATIONAL_DYADIC>(), rotmat));
93Vector<T,3> Khathatvelo (0.);
94Khathatvelo = eler::matrixVectorMultiply(Khathat, fvel - velocity);//in Schachar Beerman is K
95
96 /* / 0 1 2 \
97 * | 3 4 5 |
98 * \ 6 7 8 /
99 */
100
101Vector<T,3> F12, F13, F23, F21, F31, F32;
102Vector<T,3> force = Vector<T,3> (0.,0.,0.);
103Vector<T,9> B12, B13, B23, B21, B31, B32;
104Vector<T,9> invB12, invB13, invB23, invB21, invB31, invB32;
105B12 = Vector<T,9> (1.,0.,0.,0.,1.,0.,0.,0.,1.);
106invB12 = B12;
107B13 = Vector<T,9> (1.,0.,0.,0.,0.,-1.,0.,1.,0.);
108invB13 = Vector<T,9> (1.,0.,0.,0.,0.,1.,0.,-1.,0.);
109B21 = Vector<T,9> (0.,0.,1.,1.,0.,0.,0.,-1.,0.);
110invB21 = Vector<T,9>(0.,1.,0.,0.,0.,-1.,1.,0.,0.);
111B23 = Vector<T,9>(0.,0.,1.,1.,0.,0.,0.,1.,0.);
112invB23 = Vector<T,9>(0.,1.,0.,0.,0.,1.,1.,0.,0.);
113B31 = Vector<T,9>(0.,1.,0.,0.,0.,1.,1.,0.,0.);
114invB31 = Vector<T,9>(0.,0.,1.,1.,0.,0.,0.,1.,0.);
115B32 = Vector<T,9>(0.,0.,-1.,0.,1.,0.,1.,0.,0.);
116invB32 = Vector<T,9>(0.,0.,1.,0.,1.,0.,-1.,0.,0.);
117
118Vector<T,9> L (0.0501, 0.0329, 0.0, 0.0182, 0.0173, 0.00, 0.00, 0.00, 0.0373);
119
120
121 Vector <T,3> posVector (0.);
122 posVector = positionArray;
123 Vector<T,3> x (1.,0.,0.);
124 Vector<T,3> y (0.,1.,0.);
125 Vector<T,3> z (0.,0.,1.);
126 Vector<T,3> gradx (0.);
127 Vector<T,3> grady (0.);
128 Vector<T,3> gradz (0.);
129 veloGradientInterpFD<T,DESCRIPTOR>(gradx,blockInterpPhysVelF, x,posVector, this->_converter.getPhysDeltaX()/100.);
130 veloGradientInterpFD<T,DESCRIPTOR>(grady,blockInterpPhysVelF, y,posVector, this->_converter.getPhysDeltaX()/100.);
131 veloGradientInterpFD<T,DESCRIPTOR>(gradz,blockInterpPhysVelF, z,posVector, this->_converter.getPhysDeltaX()/100.);
132
133 /* / 0 1 2 \ vx/dx vx/dy vx/dz
134 | 3 4 5 | vy/dx vy/dy vy/dz
135 \ 6 7 8 / vz/dx vz/dy vz/dz */
136
137 Vector <T,9> gradvel (0.);
138 gradvel[0] = gradx[0];
139 gradvel[3] = gradx[1];
140 gradvel[6] = gradx[2];
141 gradvel[1] = grady[0];
142 gradvel[4] = grady[1];
143 gradvel[7] = grady[2];
144 gradvel[2] = gradz[0];
145 gradvel[5] = gradz[1];
146 gradvel[8] = gradz[2];
147
148using namespace eler;
149using namespace util;
150
151T lim = std::numeric_limits<T>::min();
152
153T Const = M_PI*M_PI*particle.template getField<PHYSPROPERTIES,RADIUS>()*particle.template getField<PHYSPROPERTIES,RADIUS>()*this->_converter.getPhysViscosity()*this->_converter.getPhysDensity()/sqrt(this->_converter.getPhysViscosity());//this->_converter.getPhysViscosity()*this->_converter.getPhysDensity()*M_PI*particle.template getField<PHYSPROPERTIES,RADIUS>();
154
155if(std::fabs(gradvel[1]) < lim )
156F12 = Vector<T,3> (0.,0.,0.);
157else
158F12 = Const*gradvel[1]/sqrt(std::fabs(gradvel[1]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B12,matrixMultiply(L,matrixMultiply(invB12,Khathat)))),fvel - velocity);
159if(std::fabs(gradvel[2]) < lim )
160F13 = Vector<T,3> (0.,0.,0.);
161else
162F13 = Const*gradvel[2]/sqrt(std::fabs(gradvel[2]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B13,matrixMultiply(L,matrixMultiply(invB13,Khathat)))),fvel - velocity) ;
163if(std::fabs(gradvel[3]) < lim )
164F21 = Vector<T,3> (0.,0.,0.);
165else
166F21 = Const*gradvel[3]/sqrt(std::fabs(gradvel[3]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B21,matrixMultiply(L,matrixMultiply(invB21,Khathat)))),fvel - velocity) ;
167if(std::fabs(gradvel[5]) < lim )
168F23 = Vector<T,3> (0.,0.,0.);
169else
170F23 = Const*gradvel[5]/sqrt(std::fabs(gradvel[5]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B23,matrixMultiply(L,matrixMultiply(invB23,Khathat)))),fvel - velocity) ;
171if(std::fabs(gradvel[6]) < lim )
172F31 = Vector<T,3> (0.,0.,0.);
173else
174F31 = Const*gradvel[6]/sqrt(std::fabs(gradvel[6]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B31,matrixMultiply(L,matrixMultiply(invB31,Khathat)))),fvel - velocity) ;
175if(std::fabs(gradvel[7]) < lim )
176F32 = Vector<T,3> (0.,0.,0.);
177else
178F32 = Const*gradvel[7]/sqrt(std::fabs(gradvel[7]))*matrixVectorMultiply(matrixMultiply(Khathat,matrixMultiply(B32,matrixMultiply(L,matrixMultiply(invB32,Khathat)))),fvel - velocity) ;
179
180
181force = F12 + F13 + F21 + F23 + F31 + F32;
182
183 particle.template setField<FORCING,FORCE> ( particle.template getField<FORCING,FORCE>()+force);
184
185 }
186}
#define M_PI
Cuboid< T, D > & getCuboid()
constexpr bool providesParallelization()
Provides group PARALLELIZATION.
std::conditional_t< DESCRIPTOR::d==2, BlockLatticeInterpPhysVelocity2D< T, DESCRIPTOR >, BlockLatticeInterpPhysVelocity3D< T, DESCRIPTOR > > BlockLatticeInterpPhysVelocity
Definition aliases.h:349
Vector(T &&t, Ts &&... ts) -> Vector< std::remove_cvref_t< T >, 1+sizeof...(Ts)>

References olb::Vector< T, Size >::data(), M_PI, olb::particles::access::providesParallelization(), and olb::Vector().

+ Here is the call graph for this function:

◆ operator()()

template<typename T , typename DESCRIPTOR , typename PARTICLETYPE >
bool olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >::operator() ( T output[],
const int input[] )
override

Definition at line 191 of file latticeSpheroidLiftForce.hh.

192{
193 using namespace descriptors;
194 using namespace particles;
195 // iterate over all particles in _indicator //TODO: add periodic treatment analogous to momentum exchange
196 forParticlesInParticleSystem<T,PARTICLETYPE,conditions::valid_particles>( _particleSystem,
197 [&](Particle<T,PARTICLETYPE>& particle){
198
199 evaluate(output, particle, particle.getId());
200
201
202 });
203
204
205 return true;
206}
void evaluate(T output[], particles::Particle< T, PARTICLETYPE > &particle, int iP)

Member Data Documentation

◆ serializeForce

template<typename T , typename DESCRIPTOR , typename PARTICLETYPE >
bool olb::BlockLatticeSpheroidLiftForce< T, DESCRIPTOR, PARTICLETYPE >::serializeForce = false
staticconstexpr

Definition at line 82 of file latticeSpheroidLiftForce.h.


The documentation for this class was generated from the following files: