29#ifndef PARTICLE_DYNAMICS_FUNCTIONS_H
30#define PARTICLE_DYNAMICS_FUNCTIONS_H
45 if constexpr (D == 3) {
63template <
typename T,
typename PARTICLECONTACTTYPE,
typename WALLCONTACTTYPE>
64void communicateContacts(ContactContainer<T, PARTICLECONTACTTYPE, WALLCONTACTTYPE>& contactContainer);
72template<
unsigned D,
typename T>
90template <
typename T,
typename PARTICLETYPE>
93 using namespace descriptors;
96 particle.template getField<GENERAL,POSITION>();
98 particle.template getField<MOBILITY,VELOCITY>();
100 angVelocity(particle.template getField<MOBILITY,ANG_VELOCITY>());
109template<
typename T,
typename PARTICLETYPE>
112 T serializedFunctorForceField[],
int iP )
114 constexpr unsigned D = PARTICLETYPE::d;
116 constexpr int serialSize = D+Drot+1;
118 const int idxForce = iP*(serialSize);
119 const int idxTorque = idxForce+D;
120 const int idxTorqueEnd = idxTorque+Drot;
124 force = std::vector<T>( serializedFunctorForceField+idxForce,
125 serializedFunctorForceField+idxTorque);
127 torque = std::vector<T>( serializedFunctorForceField+idxTorque,
128 serializedFunctorForceField+idxTorqueEnd );
133template<
typename T,
typename PARTICLETYPE>
135 T serializedFunctorForceField[],
int iP )
137 constexpr unsigned D = PARTICLETYPE::d;
138 constexpr int serialSize = D;
140 const int idxForce = iP*(serialSize);
141 const int idxForceEnd = idxForce+D;
145 force = std::vector<T>( serializedFunctorForceField+idxForce,
146 serializedFunctorForceField+idxForceEnd);
150template<
typename T,
typename PARTICLETYPE,
typename FORCEFUNCTOR>
153 constexpr unsigned D = PARTICLETYPE::d;
158 T serializedFunctorForceField[forceF.getTargetDim()];
162 if constexpr (std::is_base_of<BlockF<T,PARTICLETYPE::d>, FORCEFUNCTOR>::value){
163 for (
int iS=0; iS<forceF.getTargetDim(); ++iS) {
164 serializedFunctorForceField[iS] = 0.;
169 forceF(serializedFunctorForceField, input);
173 for (std::size_t iP=iP0; iP<particleSystem.
size(); iP++) {
174 auto particle = particleSystem.
get(iP);
177 std::size_t iPeval = iP-iP0;
178 if constexpr ( PARTICLETYPE::template providesNested<descriptors::FORCING,descriptors::TORQUE>() ) {
180 unserializeForceTorqueVoxels<T,PARTICLETYPE>( force, torque, serializedFunctorForceField, iPeval );
181 particle.template setField<descriptors::FORCING,descriptors::TORQUE>(
185 unserializeForce<T,PARTICLETYPE>( force, serializedFunctorForceField, iPeval );
191 particle.template setField<descriptors::FORCING,descriptors::FORCE>( force );
197template<
typename T,
typename PARTICLETYPE,
typename FORCEFUNCTOR,
typename PCONDITION=conditions::val
id_particles>
201 for (std::size_t iP=iP0; iP!=particleSystem.
size(); iP++) {
202 auto particle = particleSystem.
get(iP);
204 doWhenMeetingCondition<T,PARTICLETYPE,PCONDITION>( particle,
208 forceF.evaluate(output, particle, iP);
214template<
typename T,
typename PARTICLETYPE>
218 typedef std::function<bool(
const std::type_info&,
int,std::string)> FunctionType;
219 FunctionType initFunction = [](
const std::type_info& typeInfo,
int fieldSize, std::string fieldContentStr) {
224 initFunction, dynamicFieldGroups, iP );
233template<
typename T,
typename DESCRIPTOR,
typename PARTICLETYPE,
typename PARTICLECONTACTTYPE,
typename WALLCONTACTTYPE,
234 typename F=
decltype(defaults::periodicity<DESCRIPTOR::d>)>
242 F getSetupPeriodicity = defaults::periodicity<DESCRIPTOR::d>)
244 static_assert(DESCRIPTOR::template provides<descriptors::CONTACT_DETECTION>(),
245 "The field CONTACT_DETECTION must be provided.");
246 constexpr unsigned D = DESCRIPTOR::d;
247 using namespace descriptors;
255 for (std::size_t iP=0; iP<particleSystem.
size(); ++iP) {
256 auto particle = particleSystem.
get(iP);
259 particleSystem, contactContainer, iP, particle,
260 solidBoundaries, getSetupPeriodicity );
264 contact::communicateContacts<T,PARTICLECONTACTTYPE,WALLCONTACTTYPE>(contactContainer);
268template<
typename T,
typename DESCRIPTOR,
typename PARTICLETYPE,
typename PARTICLECONTACTTYPE,
typename WALLCONTACTTYPE,
269 typename F=
decltype(defaults::periodicity<DESCRIPTOR::d>)>
276 F getSetupPeriodicity = defaults::periodicity<DESCRIPTOR::d>)
278 std::vector<SolidBoundary<T,DESCRIPTOR::d>> solidBoundaries = std::vector<SolidBoundary<T,DESCRIPTOR::d>>();
284template<
typename T,
typename DESCRIPTOR,
typename PARTICLETYPE,
typename PARTICLECONTACTTYPE,
typename WALLCONTACTTYPE,
285 typename F=
decltype(defaults::periodicity<DESCRIPTOR::d>)>
293 F getSetupPeriodicity = defaults::periodicity<DESCRIPTOR::d>)
295 static_assert(DESCRIPTOR::template provides<descriptors::CONTACT_DETECTION>(),
296 "The field CONTACT_DETECTION must be provided.");
297 constexpr unsigned D = DESCRIPTOR::d;
298 using namespace descriptors;
304 communication::forParticlesInSuperParticleSystem<T, PARTICLETYPE, conditions::valid_particles>(
308 const std::size_t globalParticleID =
309 particle.template getField<PARALLELIZATION, ID>();
312 particleSystem, contactContainer,
313 globalParticleID, particle, solidBoundaries,
314 getSetupPeriodicity, globiC );
319template<
typename T,
typename DESCRIPTOR,
typename PARTICLETYPE,
typename PARTICLECONTACTTYPE,
typename WALLCONTACTTYPE,
320 typename F=
decltype(defaults::periodicity<DESCRIPTOR::d>)>
327 F getSetupPeriodicity = defaults::periodicity<DESCRIPTOR::d>)
329 std::vector<SolidBoundary<T,DESCRIPTOR::d>> solidBoundaries = std::vector<SolidBoundary<T,DESCRIPTOR::d>>();
334template<
typename T,
typename PARTICLETYPE>
337 using namespace descriptors;
338 constexpr unsigned D = PARTICLETYPE::d;
340 static_assert(D==3,
"ERROR: 2D version of calcKineticEnergy not implemented yet!");
342 T pMass = particle.template getField<PHYSPROPERTIES, MASS>();
343 Vector<T,D> vel(particle.template getField<MOBILITY, VELOCITY>());
344 Vector<T,Drot> pMofi(particle.template getField<PHYSPROPERTIES, MOFI>());
345 Vector<T,Drot> angVel(particle.template getField<MOBILITY, ANG_VELOCITY>());
347 T eTrans = T{0.5} * pMass * util::normSqr<T,D>(vel);
348 T eRot = .5*(pMofi[0]*angVel[0]*angVel[0]
349 +pMofi[1]*angVel[1]*angVel[1]
350 +pMofi[2]*angVel[2]*angVel[2]);
351 T eKin = eTrans+eRot;
Storage for dynamic field groups (Prototype for ParticleSystem)
Representation of a statistic for a parallel 2D geometry.
Super class maintaining block lattices for a cuboid decomposition.
CuboidGeometry< T, D > & getCuboidGeometry()
Read and write access to cuboid geometry.
Conversion between physical and lattice units, as well as discretization.
auto & get()
Expose container.
constexpr std::size_t size()
Size of ParticleSystem.
void communicateContacts(ContactContainer< T, PARTICLECONTACTTYPE, WALLCONTACTTYPE > &contactContainer)
void applyLocalParticleForce(FORCEFUNCTOR &forceF, ParticleSystem< T, PARTICLETYPE > &particleSystem, std::size_t iP0=0)
Apply boundary force provided by force functor to the particle center as torque and force.
void applySerializableParticleForce(FORCEFUNCTOR &forceF, ParticleSystem< T, PARTICLETYPE > &particleSystem, std::size_t iP0=0)
Apply boundary force provided by force functor to the particle center as torque and force.
void unserializeForceTorqueVoxels(Vector< T, PARTICLETYPE::d > &force, Vector< T, utilities::dimensions::convert< PARTICLETYPE::d >::rotation > &torque, T serializedFunctorForceField[], int iP)
Unserialize force field provieded by force integration functor (e.g. momentumExchange)
constexpr Vector< T, PARTICLETYPE::d > calculateLocalVelocity(Particle< T, PARTICLETYPE > &particle, const PhysR< T, PARTICLETYPE::d > &input)
T calcKineticEnergy(Particle< T, PARTICLETYPE > &particle)
void initializeParticle(DynamicFieldGroupsD< T, typename PARTICLETYPE::fields_t > &dynamicFieldGroups, std::size_t iP)
Initialize all fields in particle (necessary for clang)
void unserializeForce(Vector< T, PARTICLETYPE::d > &force, T serializedFunctorForceField[], int iP)
Unserialize force field provieded by force integration functor (e.g. stokesDragForce)
void coupleResolvedParticlesToLattice(ParticleSystem< T, PARTICLETYPE > &particleSystem, contact::ContactContainer< T, PARTICLECONTACTTYPE, WALLCONTACTTYPE > &contactContainer, const SuperGeometry< T, DESCRIPTOR::d > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, UnitConverter< T, DESCRIPTOR > const &converter, std::vector< SolidBoundary< T, DESCRIPTOR::d > > &solidBoundaries, F getSetupPeriodicity=defaults::periodicity< DESCRIPTOR::d >)
Couple particle to lattice and detect contacts of resolved particles.
constexpr bool isPeriodic(const Vector< bool, D > &periodic)
void setSuperParticleField(const SuperGeometry< T, DESCRIPTOR::d > &sGeometry, SuperLattice< T, DESCRIPTOR > &sLattice, UnitConverter< T, DESCRIPTOR > const &converter, Particle< T, PARTICLETYPE > &particle, const Vector< bool, DESCRIPTOR::d > &periodicity)
Set particle field with peridic support.
constexpr Vector< T, D > calculateLocalVelocity(const Vector< T, D > &rotationCenter, const Vector< T, D > &velocity, const Vector< T, utilities::dimensions::convert< D >::rotation > &angularVelocity, const Vector< T, D > &position)
Calculate local velocity.
Top level namespace for all of OpenLB.
auto crossProduct(const ScalarVector< T, D, IMPL > &a, const ScalarVector< T, D, IMPL_ > &b)
Traversal of nested field contents for output and initialization.
void cleanContacts()
Clean contacts - remove "empty" contacts.
static constexpr Vector< T, utilities::dimensions::convert< D >::rotation > calculate(Vector< T, D > force, PhysR< T, D > lever)
Converts dimensions by deriving from given cartesian dimension D.