24#ifndef GPU_CUDA_DYNAMICS_HH
25#define GPU_CUDA_DYNAMICS_HH
36template <
typename T,
typename DESCRIPTOR>
class DeviceContext;
37template <
typename T,
typename DESCRIPTOR>
class DataOnlyCell;
40template <
typename T,
typename DESCRIPTOR>
62 for (
unsigned iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
68 T rho, T* u, T* pi) __device__ {
70 for (
unsigned iPop=0; iPop < DESCRIPTOR::q; ++iPop) {
80template <
typename T,
typename DESCRIPTOR,
typename DYNAMICS>
87 _parameters{parameters} {
92 return DYNAMICS().apply(cell, *_parameters);
96 return DYNAMICS::MomentaF().computeRho(cell);
99 DYNAMICS::MomentaF().computeU(cell, u);
102 DYNAMICS::MomentaF().computeJ(cell, j);
105 DYNAMICS::MomentaF().computeRhoU(cell, rho, u);
108 DYNAMICS::MomentaF().computeStress(cell, rho, u, pi);
111 DYNAMICS::MomentaF().computeAllMomenta(cell, rho, u, pi);
115 if constexpr (DYNAMICS::parameters::template contains<descriptors::OMEGA>()) {
116 return _parameters->template get<descriptors::OMEGA>();
123 return DYNAMICS().computeEquilibrium(iPop, rho, u);
127 typename DYNAMICS::MomentaF().defineRho(cell, rho);
131 typename DYNAMICS::MomentaF().defineU(cell, u);
135 typename DYNAMICS::MomentaF().defineRhoU(cell, rho, u);
139 typename DYNAMICS::MomentaF().defineAllMomenta(cell, rho, u, pi);
143 typename DYNAMICS::MomentaF().inverseShiftRhoU(cell, rho, u);
149 template <
typename T,
typename DESCRIPTOR>
152 collisionO->collide(lattice, iCell);
158 template <
typename T,
typename DESCRIPTOR>
161 statistic = collisionO->collide(lattice, iCell);
Implementation of gpu::cuda::Dynamics for concrete DYNAMICS.
void defineU(DataOnlyCell< T, DESCRIPTOR > &cell, T *u) override __device__
CellStatistic< T > collide(DeviceContext< T, DESCRIPTOR > lattice, CellID iCell) override __device__
void computeAllMomenta(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) override __device__
void computeJ(DataOnlyCell< T, DESCRIPTOR > &cell, T *j) override __device__
void inverseShiftRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) override __device__
T computeEquilibrium(int iPop, T rho, T *u) override __device__
ConcreteDynamics(ParametersOfOperatorD< T, DESCRIPTOR, DYNAMICS > *parameters) __device__
void defineAllMomenta(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) override __device__
void computeU(DataOnlyCell< T, DESCRIPTOR > &cell, T *u) override __device__
void defineRho(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho) override __device__
void computeStress(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) override __device__
T computeRho(DataOnlyCell< T, DESCRIPTOR > &cell) override __device__
T getOmegaOrFallback(T fallback) override __device__
void defineRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) override __device__
void computeRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) override __device__
Device-side implementation of the data-only Cell concept for collision steps.
Structure for passing pointers to on-device data into CUDA kernels.
Top level namespace for all of OpenLB.
std::uint32_t CellID
Type for sequential block-local cell indices.
typename ParametersD< T, DESCRIPTOR >::template include< typename OPERATOR::parameters > ParametersOfOperatorD
Deduce ParametersD of OPERATOR w.r.t. T and DESCRIPTOR.
Return value of any collision.
On-device field mirroring BlockDynamicsMap.
Last node in a MaskedDynamics chain in kernel::call_operators.
bool operator()(DeviceContext< T, DESCRIPTOR > &lattice, CellID iCell, CellStatistic< T > &statistic) __device__
bool operator()(DeviceContext< T, DESCRIPTOR > &lattice, CellID iCell) __device__
Virtual interface for device-side dynamically-dispatched dynamics access.
void iniRegularized(DataOnlyCell< T, DESCRIPTOR > &cell, T rho, T *u, T *pi) __device__
virtual T getOmegaOrFallback(T fallback) __device__=0
virtual void defineRho(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho) __device__=0
virtual T computeRho(DataOnlyCell< T, DESCRIPTOR > &cell) __device__=0
virtual void computeStress(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) __device__=0
virtual void defineAllMomenta(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) __device__=0
virtual void computeRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) __device__=0
virtual void inverseShiftRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) __device__=0
virtual void computeAllMomenta(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u, T *pi) __device__=0
virtual T computeEquilibrium(int iPop, T rho, T *u) __device__=0
virtual void computeJ(DataOnlyCell< T, DESCRIPTOR > &cell, T *j) __device__=0
virtual CellStatistic< T > collide(DeviceContext< T, DESCRIPTOR > lattice, CellID iCell) __device__=0
virtual void defineRhoU(DataOnlyCell< T, DESCRIPTOR > &cell, T &rho, T *u) __device__=0
virtual void computeU(DataOnlyCell< T, DESCRIPTOR > &cell, T *u) __device__=0
void iniEquilibrium(DataOnlyCell< T, DESCRIPTOR > &cell, T rho, T *u) __device__
virtual void defineU(DataOnlyCell< T, DESCRIPTOR > &cell, T *u) __device__=0