OpenLB 1.7
Loading...
Searching...
No Matches
operator.h
Go to the documentation of this file.
1/* This file is part of the OpenLB library
2 *
3 * Copyright (C) 2022 Adrian Kummerlaender
4 * E-mail contact: info@openlb.net
5 * The most recent release of OpenLB can be downloaded at
6 * <http://www.openlb.net/>
7 *
8 * This program is free software; you can redistribute it and/or
9 * modify it under the terms of the GNU General Public License
10 * as published by the Free Software Foundation; either version 2
11 * of the License, or (at your option) any later version.
12 *
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public
19 * License along with this program; if not, write to the Free
20 * Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
21 * Boston, MA 02110-1301, USA.
22*/
23
24#ifndef GPU_CUDA_OPERATOR_H
25#define GPU_CUDA_OPERATOR_H
26
27#include "core/operator.h"
28
29#include "dynamics/dynamics.h"
30
31#include "device.h"
32#include "mask.h"
33#include "context.h"
34#include "dynamics.h"
35
36namespace olb {
37
38
40template <typename T, typename DESCRIPTOR, typename DYNAMICS>
41class ConcreteBlockCollisionO<T,DESCRIPTOR,Platform::GPU_CUDA,DYNAMICS> final
42 : public BlockCollisionO<T,DESCRIPTOR,Platform::GPU_CUDA> {
43private:
44 std::unique_ptr<DYNAMICS> _dynamics;
45
48
49 gpu::cuda::Dynamics<T,DESCRIPTOR>** _dynamicsOfCells;
50
52
54 bool _modified;
55
57
64
65public:
67
68 std::type_index id() const override
69 {
70 return typeid(DYNAMICS);
71 }
72
73 std::size_t weight() const override
74 {
75 return _mask->weight();
76 }
77
78 void set(CellID iCell, bool state, bool overlap) override
79 {
81 if constexpr (!std::is_same_v<DYNAMICS,NoDynamics<T,DESCRIPTOR>>) {
82 if (!overlap) {
83 _mask->set(iCell, state);
84 }
85 }
86 if (state) {
87 _dynamicsOfCells[iCell] = _deviceDynamics.get();
88 }
89 _modified = true;
90 }
91
93 {
94 return _dynamics.get();
95 }
96
98
100
105 CollisionDispatchStrategy strategy) override
106 {
107 switch (strategy) {
109 return applyDominant(block, subdomain);
111 return applyIndividual(block, subdomain);
112 default:
113 throw std::runtime_error("Invalid collision dispatch strategy");
114 }
115 }
116
117};
118
119
121template <typename T, typename DESCRIPTOR, CONCEPT(CellOperator) OPERATOR>
122class ConcreteBlockO<T,DESCRIPTOR,Platform::GPU_CUDA,OPERATOR,OperatorScope::PerCell> final
123 : public BlockO<T,DESCRIPTOR,Platform::GPU_CUDA> {
124private:
126 bool _modified;
127
129
130public:
132
133 std::type_index id() const override
134 {
135 return typeid(OPERATOR);
136 }
137
138 void set(CellID iCell, bool state) override
139 {
140 if (state) {
141 _cells.push_back(iCell);
142 _modified = true;
143 }
144 }
145
147 {
148 _modified = false;
149 }
150
152
153};
154
155
157template <typename T, typename DESCRIPTOR, CONCEPT(CellOperator) OPERATOR>
159 : public BlockO<T,DESCRIPTOR,Platform::GPU_CUDA> {
160private:
162 bool _modified;
163
165
167
168public:
170
171 std::type_index id() const override
172 {
173 return typeid(OPERATOR);
174 }
175
176 void set(CellID iCell, bool state) override
177 {
178 if (state) {
179 _cells.push_back(iCell);
180 _modified = true;
181 }
182 }
183
185 {
186 _modified = false;
187 _parameters = &block.template getData<OperatorParameters<OPERATOR>>();
188 }
189
191
192};
193
194
196
199template <typename T, typename DESCRIPTOR, CONCEPT(BlockOperator) OPERATOR>
200class ConcreteBlockO<T,DESCRIPTOR,Platform::GPU_CUDA,OPERATOR,OperatorScope::PerBlock> final
201 : public BlockO<T,DESCRIPTOR,Platform::GPU_CUDA> {
202public:
203 std::type_index id() const override
204 {
205 return typeid(OPERATOR);
206 }
207
208 void set(CellID iCell, bool state) override
209 {
210 throw std::logic_error("BlockO::set not supported for OperatorScope::PerBlock");
211 }
212
214 {
215 OPERATOR().setup(block);
216 }
217
219 {
220 OPERATOR().apply(block);
221 }
222
223};
224
225
227template <typename COUPLER, typename COUPLEES>
229 final : public AbstractCouplingO<COUPLEES> {
230private:
231 template <typename VALUED_DESCRIPTOR>
232 using ptr_to_lattice = ConcreteBlockLattice<typename VALUED_DESCRIPTOR::value_t,
233 typename VALUED_DESCRIPTOR::descriptor_t,
235
236 utilities::TypeIndexedTuple<typename COUPLEES::template map_values<
238 >> _lattices;
239
240 typename AbstractCouplingO<COUPLEES>::ParametersD _parameters;
241
242 std::unique_ptr<ConcreteBlockMask<typename COUPLEES::values_t::template get<0>::value_t,
243 Platform::GPU_CUDA>> _mask;
244
245public:
246 template <typename LATTICES>
247 ConcreteBlockCouplingO(LATTICES&& lattices):
248 _lattices{lattices}
249 { }
250
251 std::type_index id() const override {
252 return typeid(COUPLER);
253 }
254
256 return _parameters;
257 }
258
259 void set(CellID iCell, bool state) override;
260
261 void execute() override;
262
263};
264
266template <typename COUPLER, typename COUPLEES>
268 final : public AbstractCouplingO<COUPLEES> {
269private:
270 template <typename VALUED_DESCRIPTOR>
271 using ptr_to_lattice = ConcreteBlockLattice<typename VALUED_DESCRIPTOR::value_t,
272 typename VALUED_DESCRIPTOR::descriptor_t,
274
275 utilities::TypeIndexedTuple<typename COUPLEES::template map_values<
277 >> _lattices;
278
279 typename COUPLER::parameters::template decompose_into<
281 > _parameters;
282
283 std::unique_ptr<ConcreteBlockMask<typename COUPLEES::values_t::template get<0>::value_t,
284 Platform::GPU_CUDA>> _mask;
285
286public:
287 template <typename LATTICES>
288 ConcreteBlockCouplingO(LATTICES&& lattices):
289 _lattices{lattices}
290 { }
291
292 std::type_index id() const override {
293 return typeid(COUPLER);
294 }
295
297 return _parameters;
298 }
299
300 void set(CellID iCell, bool state) override;
301
302 void execute() override;
303
304};
305
306}
307
308#include "communicator.h"
309#include "statistics.h"
310
311#endif
std::size_t weight() const override
Returns number of assigned cells.
Definition operator.h:73
void set(CellID iCell, bool state, bool overlap) override
Set whether iCell is covered by the present collision step.
Definition operator.h:78
void apply(ConcreteBlockLattice< T, DESCRIPTOR, Platform::GPU_CUDA > &block, ConcreteBlockMask< T, Platform::GPU_CUDA > &subdomain, CollisionDispatchStrategy strategy) override
Apply collision to subdomain of block using strategy.
Definition operator.h:103
Collision operation of concrete DYNAMICS on concrete block lattices of PLATFORM.
AbstractCouplingO< COUPLEES >::AbstractParameters & getParameters() override
Return reference to parameters of coupling operator.
Definition operator.h:296
AbstractCouplingO< COUPLEES >::AbstractParameters & getParameters() override
Return reference to parameters of coupling operator.
Definition operator.h:255
Coupling of COUPLEES using concrete OPERATOR with SCOPE on PLATFORM lattices.
Definition operator.h:142
Implementation of BlockLattice on a concrete PLATFORM.
void setup(ConcreteBlockLattice< T, DESCRIPTOR, Platform::GPU_CUDA > &block) override
Definition operator.h:184
void set(CellID iCell, bool state) override
Set whether iCell is covered by the operator (optional)
Definition operator.h:176
void set(CellID iCell, bool state) override
Set whether iCell is covered by the operator (optional)
Definition operator.h:208
void apply(ConcreteBlockLattice< T, DESCRIPTOR, Platform::GPU_CUDA > &block) override
Definition operator.h:218
void setup(ConcreteBlockLattice< T, DESCRIPTOR, Platform::GPU_CUDA > &block) override
Definition operator.h:213
void setup(ConcreteBlockLattice< T, DESCRIPTOR, Platform::GPU_CUDA > &block) override
Definition operator.h:146
void set(CellID iCell, bool state) override
Set whether iCell is covered by the operator (optional)
Definition operator.h:138
Block application of concrete OPERATOR called using SCOPE on PLATFORM.
Definition operator.h:65
Plain column for CUDA GPU targets.
Definition column.h:49
void push_back(T value)
Definition column.hh:102
Basic wrapper for device stream.
Definition device.h:121
Managed pointer for device-side memory.
Definition device.h:79
Top level namespace for all of OpenLB.
std::uint32_t CellID
Type for sequential block-local cell indices.
Platform
OpenLB execution targets.
Definition platform.h:36
@ GPU_CUDA
Vector CPU (AVX2 / AVX-512 collision)
CollisionDispatchStrategy
Collision dispatch strategy.
Definition operator.h:88
@ Individual
Apply all dynamics individually (async for Platform::GPU_CUDA)
@ Dominant
Apply dominant dynamics using mask and fallback to virtual dispatch for others.
@ 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...
Base of block-wide coupling operators executed by SuperLatticeCoupling.
Definition operator.h:115
Dynamic access interface for FIELD-valued parameters.
Collision operation on concrete blocks of PLATFORM.
Definition operator.h:97
Base of block-wide operators such as post processors.
Definition operator.h:41
Concrete storage of ParametersD in olb::Data.
Interface for per-cell dynamics.
Definition interface.h:54
Set of FIELD-valued parameters.
Virtual interface for device-side dynamically-dispatched dynamics access.
Definition dynamics.hh:41
Mapping between KEYs and instances of type VALUEs.