24#ifndef SISD_OPERATOR_H
25#define SISD_OPERATOR_H
45template <
typename T,
typename DESCRIPTOR,
typename DYNAMICS>
52 _parameters{parameters} {
56 return DYNAMICS().apply(cell, *_parameters);
60 return typename DYNAMICS::MomentaF().computeRho(cell);
63 typename DYNAMICS::MomentaF().computeU(cell, u);
66 typename DYNAMICS::MomentaF().computeJ(cell, j);
69 typename DYNAMICS::MomentaF().computeRhoU(cell, rho, u);
72 typename DYNAMICS::MomentaF().computeStress(cell, rho, u, pi);
75 typename DYNAMICS::MomentaF().computeAllMomenta(cell, rho, u, pi);
79 if constexpr (DYNAMICS::parameters::template contains<descriptors::OMEGA>()) {
80 return _parameters->template get<descriptors::OMEGA>();
84 __builtin_unreachable();
88 return DYNAMICS().computeEquilibrium(iPop, rho, u);
92 typename DYNAMICS::MomentaF().defineRho(cell, rho);
96 typename DYNAMICS::MomentaF().defineU(cell, u);
100 typename DYNAMICS::MomentaF().defineRhoU(cell, rho, u);
104 typename DYNAMICS::MomentaF().defineAllMomenta(cell, rho, u, pi);
108 typename DYNAMICS::MomentaF().inverseShiftRhoU(cell, rho, u);
124template <
typename T,
typename DESCRIPTOR,
typename DYNAMICS>
128 std::unique_ptr<DYNAMICS> _dynamics;
129 std::unique_ptr<cpu::Dynamics<T,DESCRIPTOR,Platform::CPU_SISD>> _concreteDynamics;
136 std::vector<CellID> _cells;
146 auto& parameters = *_parameters;
149 #ifdef PARALLEL_MODE_OMP
150 #pragma omp declare reduction(+ : typename LatticeStatistics<T>::Aggregatable : omp_out += omp_in) initializer (omp_priv={})
153 if constexpr (DESCRIPTOR::d == 3) {
154 #ifdef PARALLEL_MODE_OMP
155 #pragma omp parallel for schedule(dynamic,1) reduction(+ : statistics)
157 for (
int iX=0; iX < block.
getNx(); ++iX) {
158 for (
int iY=0; iY < block.
getNy(); ++iY) {
159 std::size_t iCell = block.
getCellId(iX,iY,0);
160 for (
int iZ=0; iZ < block.
getNz(); ++iZ) {
162 if (
auto cellStatistic = mask[iCell] ? DYNAMICS().apply(cell, parameters)
163 : _dynamicsOfCells[iCell]->
collide(cell)) {
164 statistics.increment(cellStatistic.rho, cellStatistic.uSqr);
171 #ifdef PARALLEL_MODE_OMP
172 #pragma omp parallel for schedule(dynamic,1) reduction(+ : statistics)
174 for (
int iX=0; iX < block.
getNx(); ++iX) {
175 std::size_t iCell = block.
getCellId(iX,0);
176 for (
int iY=0; iY < block.
getNy(); ++iY) {
178 if (
auto cellStatistic = mask[iCell] ? DYNAMICS().apply(cell, parameters)
179 : _dynamicsOfCells[iCell]->
collide(cell)) {
180 statistics.increment(cellStatistic.rho, cellStatistic.uSqr);
198 if (_mask->operator[](iCell)) {
199 _cells.push_back(iCell);
205 auto& parameters = *_parameters;
208 #ifdef PARALLEL_MODE_OMP
209 #pragma omp declare reduction(+ : typename LatticeStatistics<T>::Aggregatable : omp_out += omp_in) initializer (omp_priv={})
212 #ifdef PARALLEL_MODE_OMP
213 #pragma omp parallel for schedule(static) reduction(+ : statistics)
215 for (std::size_t i=0; i < _cells.size(); ++i) {
216 std::size_t iCell = _cells[i];
218 if (
auto cellStatistic = DYNAMICS().apply(cell, parameters)) {
219 statistics.increment(cellStatistic.rho, cellStatistic.uSqr);
228 _dynamics(new DYNAMICS()),
229 _parameters(nullptr),
235 std::type_index
id()
const override
237 return typeid(DYNAMICS);
242 return _mask->weight();
245 void set(
CellID iCell,
bool state,
bool overlap)
override
248 if constexpr (!std::is_same_v<DYNAMICS,NoDynamics<T,DESCRIPTOR>>) {
250 _mask->set(iCell, state);
254 _dynamicsOfCells[iCell] = _concreteDynamics.get();
261 return _dynamics.get();
266 _parameters = &block.template getData<OperatorParameters<DYNAMICS>>().parameters;
267 _mask = &block.template getData<DynamicsMask<DYNAMICS>>();
268 if constexpr (dynamics::has_parametrized_momenta_v<DYNAMICS>) {
269 _dynamics->setMomentaParameters(_parameters);
274 _dynamicsOfCells = block.template getField<cpu::DYNAMICS<T,DESCRIPTOR,Platform::CPU_SISD>>()[0].data();
287 return applyDominant(block, subdomain);
289 return applyIndividual(block, subdomain);
291 throw std::runtime_error(
"Invalid collision dispatch strategy");
299template <
typename T,
typename DESCRIPTOR, CONCEPT(CellOperator) OPERATOR>
301 :
public BlockO<T,DESCRIPTOR,Platform::CPU_SISD> {
303 std::vector<CellID> _cells;
307 std::type_index
id()
const override
309 return typeid(OPERATOR);
315 _cells.emplace_back(iCell);
326 std::sort(_cells.begin(), _cells.end());
327 _cells.erase(std::unique(_cells.begin(), _cells.end()), _cells.end());
330 if (_cells.size() > 0) {
332 #ifdef PARALLEL_MODE_OMP
333 #pragma omp parallel for schedule(static) firstprivate(cell)
335 for (
CellID iCell : _cells) {
337 OPERATOR().apply(cell);
345template <
typename T,
typename DESCRIPTOR, CONCEPT(CellOperator) OPERATOR>
347 :
public BlockO<T,DESCRIPTOR,Platform::CPU_SISD> {
349 std::vector<CellID> _cells;
355 std::type_index
id()
const override
357 return typeid(OPERATOR);
363 _cells.emplace_back(iCell);
370 _parameters = &block.template getData<OperatorParameters<OPERATOR>>().parameters;
376 std::sort(_cells.begin(), _cells.end());
377 _cells.erase(std::unique(_cells.begin(), _cells.end()), _cells.end());
380 if (_cells.size() > 0) {
382 #ifdef PARALLEL_MODE_OMP
383 #pragma omp parallel for schedule(static) firstprivate(cell)
385 for (
CellID iCell : _cells) {
387 OPERATOR().apply(cell, *_parameters);
399template <
typename T,
typename DESCRIPTOR, CONCEPT(BlockOperator) OPERATOR>
401 :
public BlockO<T,DESCRIPTOR,Platform::CPU_SISD> {
403 std::type_index
id()
const override
405 return typeid(OPERATOR);
410 throw std::logic_error(
"BlockO::set not supported for OperatorScope::PerBlock");
415 OPERATOR().setup(block);
420 OPERATOR().apply(block);
LatticeStatistics< T > & getStatistics()
Return a handle to the LatticeStatistics object.
std::size_t getNcells() const
Get number of cells.
int getNy() const
Read only access to block height.
int getNx() const
Read only access to block width.
CellID getCellId(LatticeR< D > latticeR) const
Get 1D cell ID.
int getNz() const
Read only access to block height.
Collision operation of concrete DYNAMICS on concrete block lattices of PLATFORM.
Implementation of BlockLattice on a concrete PLATFORM.
Block application of concrete OPERATOR called using SCOPE on PLATFORM.
Cell concept for concrete block lattices on CPU platforms.
void setCellId(std::size_t iCell)
Implementation of cpu::Dynamics for concrete DYNAMICS on SISD blocks.
void defineU(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T *u) override
void defineAllMomenta(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u, T *pi) override
void computeStress(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u, T *pi) override
T computeRho(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell) override
ConcreteDynamics(ParametersOfOperatorD< T, DESCRIPTOR, DYNAMICS > *parameters)
CellStatistic< T > collide(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell) override
void inverseShiftRhoU(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u) override
T computeEquilibrium(int iPop, T rho, T *u) override
void defineRhoU(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u) override
T getOmegaOrFallback(T fallback) override
void computeU(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T *u) override
void computeRhoU(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u) override
void defineRho(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho) override
void computeJ(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T *j) override
void computeAllMomenta(cpu::Cell< T, DESCRIPTOR, Platform::CPU_SISD > &cell, T &rho, T *u, T *pi) override
Interface for post-processing steps – header file.
Top level namespace for all of OpenLB.
std::uint32_t CellID
Type for sequential block-local cell indices.
Platform
OpenLB execution targets.
CollisionDispatchStrategy
Collision dispatch strategy.
@ Individual
Apply all dynamics individually (async for Platform::GPU_CUDA)
@ Dominant
Apply dominant dynamics using mask and fallback to virtual dispatch for others.
typename ParametersD< T, DESCRIPTOR >::template include< typename OPERATOR::parameters > ParametersOfOperatorD
Deduce ParametersD of OPERATOR w.r.t. T and DESCRIPTOR.
@ PerBlock
Per-block application, i.e. OPERATOR::apply is passed a ConcreteBlockLattice.
@ PerCell
Per-cell application, i.e. OPERATOR::apply is passed a CELL concept implementation.
@ PerCellWithParameters
Per-cell application with parameters, i.e. OPERATOR::apply is passed a CELL concept implementation an...
Collision operation on concrete blocks of PLATFORM.
Base of block-wide operators such as post processors.
Return value of any collision.
Interface for per-cell dynamics.
CPU specific field mirroring BlockDynamicsMap.
Virtual interface for dynamically-dispatched dynamics access on CPU targets.
virtual CellStatistic< T > collide(Cell< T, DESCRIPTOR, PLATFORM > &cell)=0