OpenLB 1.7
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Attributes | List of all members
olb::SuperLattice< T, DESCRIPTOR > Class Template Reference

Super class maintaining block lattices for a cuboid decomposition. More...

#include <superLattice.h>

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

Public Types

using value_t = T
 Base value type of the lattice.
 
using descriptor_t = DESCRIPTOR
 Descriptor / discrete velocity set of the lattice.
 
using block_t = ConcretizableBlockLattice<T,DESCRIPTOR>
 
- Public Types inherited from olb::SuperStructure< T, DESCRIPTOR::d >
using value_t
 

Public Member Functions

 SuperLattice (SuperGeometry< T, DESCRIPTOR::d > &superGeometry)
 Construct lattice for the cuboid decomposition of superGeometry.
 
 SuperLattice (const SuperLattice &)=delete
 
 ~SuperLattice ()=default
 
BlockLattice< T, DESCRIPTOR > & getBlock (int locC)
 Return BlockLattice with local index locC.
 
template<typename BLOCK >
BLOCK & getBlock (int locC)
 Return locC-th block lattice casted to BLOCK.
 
const BlockLattice< T, DESCRIPTOR > & getBlock (int locC) const
 Return read-only BlockLattice with local index locC.
 
template<typename BLOCK >
const BLOCK & getBlock (int locIC) const
 Return locC-th block lattice casted to BLOCK (read-only)
 
template<Platform PLATFORM, typename F >
void forBlocksOnPlatform (F f)
 Apply f to every ConcreteBlockLattice of PLATFORM.
 
void setProcessingContext (ProcessingContext context)
 Set processing context of block lattices.
 
template<typename FIELD_TYPE >
void setProcessingContext (ProcessingContext context)
 Set processing context of FIELD_TYPE in block lattices.
 
template<typename STAGE >
SuperCommunicator< T, SuperLattice > & getCommunicator (STAGE stage=STAGE())
 Return communicator for given communication stage.
 
void communicate () override
 Perform full overlap communication if needed.
 
LatticeStatistics< T > & getStatistics ()
 Return a handle to the LatticeStatistics object.
 
LatticeStatistics< T > const & getStatistics () const
 Return a constant handle to the LatticeStatistics object.
 
Cell< T, DESCRIPTOR > get (LatticeR< DESCRIPTOR::d+1 > latticeR)
 Get local cell interface.
 
template<typename... R>
std::enable_if_t< sizeof...(R)==DESCRIPTOR::d+1, Cell< T, DESCRIPTOR > > get (R... latticeR)
 Get local cell interface.
 
void initialize ()
 Initialize lattice to be ready for simulation.
 
template<typename DYNAMICS >
void defineDynamics (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator)
 Set dynamics of indicated cells to DYNAMICS.
 
template<typename DYNAMICS >
void defineDynamics (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material)
 Set dynamics of material cells to DYNAMICS.
 
template<template< typename... > typename DYNAMICS>
void defineDynamics (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator)
 Set dynamics of indicated cells to DYNAMICS<T,DESCRIPTOR>
 
template<template< typename... > typename DYNAMICS>
void defineDynamics (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material)
 Set dynamics of indicated cells to DYNAMICS<T,DESCRIPTOR>
 
void defineDynamics (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, Dynamics< T, DESCRIPTOR > *dynamics)
 Define the dynamics on a domain described by an indicator (legacy)
 
void defineDynamics (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, Dynamics< T, DESCRIPTOR > *dynamics)
 Define the dynamics on a domain with a particular material number (legacy)
 
void defineRho (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&, AnalyticalF< DESCRIPTOR::d, T, T > &rho)
 Define rho on a domain described by an indicator.
 
void defineRho (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &rho)
 Define rho on a domain with a particular material number.
 
void defineU (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Define u on a domain described by an indicator.
 
void defineU (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Define u on a domain with a particular material number.
 
void defineRhoU (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Define rho and u on a domain described by an indicator.
 
void defineRhoU (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Define rho and u on a domain with a particular material number.
 
void definePopulations (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &Pop)
 Define a population on a domain described by an indicator.
 
void definePopulations (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &Pop)
 Define a population on a domain with a particular material number.
 
void definePopulations (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, SuperF< DESCRIPTOR::d, T, T > &Pop)
 
void definePopulations (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, SuperF< DESCRIPTOR::d, T, T > &Pop)
 Define a population on a domain with a particular material number.
 
template<typename FIELD >
void defineField (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, FunctorPtr< SuperF< DESCRIPTOR::d, T, T > > &&field)
 Define an external field on a domain described by an indicator.
 
template<typename FIELD >
void defineField (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &field)
 Defines a field on a domain described by an indicator.
 
template<typename FIELD >
void defineField (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, FunctorPtr< SuperF< DESCRIPTOR::d, T, T > > &&field)
 Define an external field on a domain with a particular material number.
 
template<typename FIELD >
void defineField (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &field)
 Defines a field on a domain with a particular material number.
 
template<typename FIELD >
void defineField (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, IndicatorF2D< T > &indicator, AnalyticalF< DESCRIPTOR::d, T, T > &field)
 Defines a field on a indicated domain.
 
template<typename FIELD >
void defineField (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, IndicatorF3D< T > &indicator, AnalyticalF< DESCRIPTOR::d, T, T > &field)
 
template<typename PARAMETER >
void setParameter (FieldD< T, DESCRIPTOR, PARAMETER > field)
 Update PARAMETER in all dynamics and post processors.
 
template<typename PARAMETER , typename DYNAMICS >
void setParameterOfDynamics (FieldD< T, DESCRIPTOR, PARAMETER > &&field)
 Update PARAMETER in DYNAMICS.
 
template<typename PARAMETER , template< typename... > typename DYNAMICS>
void setParameterOfDynamics (FieldD< T, DESCRIPTOR, PARAMETER > &&field)
 Update PARAMETER in DYNAMICS<T,DESCRIPTOR>
 
template<typename PARAMETER , Platform PLATFORM, typename FIELD >
void setParameter (FieldArrayD< T, DESCRIPTOR, PLATFORM, FIELD > &fieldArray)
 Set PARAMETER to column pointers of given field array.
 
void iniEquilibrium (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Initialize by equilibrium on a domain described by an indicator.
 
void iniEquilibrium (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
 Initialize by equilibrium on a domain with a particular material number.
 
void iniEquilibrium (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, SuperF< DESCRIPTOR::d, T, T > &u)
 
void iniEquilibrium (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &rho, SuperF< DESCRIPTOR::d, T, T > &u)
 
void iniRegularized (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u, AnalyticalF< DESCRIPTOR::d, T, T > &pi)
 Initialize by non- and equilibrium on a domain described by an indicator.
 
void iniRegularized (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u, AnalyticalF< DESCRIPTOR::d, T, T > &pi)
 Initialize by non- and equilibrium on a domain with a particular material number.
 
void collideAndStream ()
 Core implementation of a single iteration of the collide and stream loop.
 
void stripeOffDensityOffset (T offset)
 Subtract constant offset from the density.
 
void statisticsOn ()
 Switch Statistics on (default on)
 
void statisticsOff ()
 Switch Statistics off (default on)
 
template<typename STAGE = stage::PostStream>
void addPostProcessor (PostProcessorGenerator< T, DESCRIPTOR > const &ppGen)
 Add a non-local post-processing step (legacy)
 
template<typename STAGE = stage::PostStream>
void addPostProcessor (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, PostProcessorGenerator< T, DESCRIPTOR > const &ppGen)
 Add a non-local post-processing step (legacy)
 
template<typename STAGE = stage::PostStream>
void addPostProcessor (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, PostProcessorGenerator< T, DESCRIPTOR > const &ppGen)
 Add a non-local post-processing step (legacy)
 
template<typename STAGE = stage::PostStream>
void addPostProcessor (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, PostProcessorPromise< T, DESCRIPTOR > &&promise)
 Add a non-local post-processing step (legacy)
 
template<typename STAGE = stage::PostStream>
void addPostProcessor (PostProcessorPromise< T, DESCRIPTOR > &&promise)
 Add a non-local post-processing step (legacy)
 
template<typename STAGE = stage::PostStream>
void executePostProcessors (STAGE stage=STAGE())
 Executes post processors for STAGE.
 
template<typename PARTNER_DESCRIPTOR >
void addLatticeCoupling (LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices)
 Adds a coupling generator for a vector of partner superLattice (legacy)
 
template<typename PARTNER_DESCRIPTOR >
void addLatticeCoupling (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices)
 Adds a coupling generator for a vector of partner superLattice (legacy)
 
template<typename PARTNER_DESCRIPTOR >
void addLatticeCoupling (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices)
 Adds a coupling generator for a vector of partner superLattice (legacy)
 
template<typename... PARTNER_DESCRIPTORS>
void addLatticeCoupling (LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, PARTNER_DESCRIPTORS &... partnerLattices)
 Adds a coupling generator for a multiple partner superLattice (legacy)
 
template<typename... PARTNER_DESCRIPTORS>
void addLatticeCoupling (FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, PARTNER_DESCRIPTORS &... partnerLattices)
 Adds a coupling generator for a multiple partner superLattice (legacy)
 
template<typename... PARTNER_DESCRIPTORS>
void addLatticeCoupling (SuperGeometry< T, DESCRIPTOR::d > &sGeometry, int material, LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, PARTNER_DESCRIPTORS &... partnerLattices)
 Adds a coupling generator for a multiple partner superLattice (legacy)
 
void executeCoupling ()
 Executes coupling with partner lattices (legacy)
 
template<typename STAGE >
void addCustomTask (std::function< void()> f)
 Schedules f for execution during every invocation of STAGE.
 
template<typename STAGE >
void executeCustomTasks (STAGE stage=STAGE())
 Executes custom tasks assigned to STAGE.
 
template<typename STAGE , typename F >
void scheduleBackgroundTask (F &&f)
 Schedule F for one-off background execution to be finished by STAGE.
 
template<typename F >
void scheduleBackgroundOutput (F &&f)
 Schedule F(iC) for one-off background output of block data.
 
template<typename CONTEXT >
void scheduleBackgroundOutputVTK (CONTEXT &&vtkContext)
 Schedule one-off background output of given VTK CONTEXT.
 
template<typename STAGE >
void waitForBackgroundTasks (STAGE stage=STAGE{})
 Block until all background tasks scheduled for STAGE are completed.
 
std::size_t getNblock () const override
 Number of data blocks for the serializable interface.
 
std::size_t getSerializableSize () const override
 Binary size for the serializer.
 
bool * getBlock (std::size_t iBlock, std::size_t &sizeBlock, bool loadingMode) override
 Return a pointer to the memory of the current block and its size for the serializable interface.
 
void postLoad () override
 
template<typename STAGE >
SuperCommunicator< T, SuperLattice< T, DESCRIPTOR > > & getCommunicator (STAGE stage)
 
- Public Member Functions inherited from olb::SuperStructure< T, DESCRIPTOR::d >
virtual ~SuperStructure ()
 Virtual Destructor for inheritance.
 
 SuperStructure (CuboidGeometry< T, D > &cuboidGeometry, LoadBalancer< T > &loadBalancer, int overlap=2)
 Construction of a super structure.
 
 SuperStructure (int overlap=1)
 Default Constructor for empty SuperStructure.
 
CuboidGeometry< T, D > & getCuboidGeometry ()
 Read and write access to cuboid geometry.
 
CuboidGeometry< T, D > const & getCuboidGeometry () const
 Read only access to cuboid geometry.
 
int getOverlap ()
 Read and write access to the overlap.
 
int getOverlap () const
 Read only access to the overlap.
 
LoadBalancer< T > & getLoadBalancer ()
 Read and write access to the load balancer.
 
LoadBalancer< T > const & getLoadBalancer () const
 Read only access to the load balancer.
 
void forCorePhysLocations (F f) const
 Iterate over discrete physical locations.
 
void forCorePhysLocations (PhysR< T, D > min, PhysR< T, D > max, F f) const
 Iterate over discrete physical locations between min and max.
 
void forCoreSpatialLocations (F f) const
 Iterate over spatial locations NOTE: Based on physical locations (as opposed to its blockStructure version)
 
void forCoreSpatialLocations (PhysR< T, D > min, PhysR< T, D > max, F f) const
 Iterate over spatial locations between min and max NOTE: Based on physical locations (as opposed to its blockStructure version)
 
- Public Member Functions inherited from olb::Serializable
virtual ~Serializable ()=default
 
template<bool includeLogOutputDir = true>
bool save (std::string fileName="", const bool enforceUint=false)
 Save Serializable into file fileName
 
template<bool includeLogOutputDir = true>
bool load (std::string fileName="", const bool enforceUint=false)
 Load Serializable from file fileName
 
bool save (std::uint8_t *buffer)
 Save Serializable into buffer of length getSerializableSize
 
bool load (const std::uint8_t *buffer)
 Load Serializable from buffer of length getSerializableSize
 

Static Public Attributes

static constexpr unsigned d = DESCRIPTOR::d
 

Additional Inherited Members

- Protected Member Functions inherited from olb::BufferSerializable
template<typename DataType >
void registerSerializable (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, DataType &data, const bool loadingMode=false)
 Register Serializable object of dynamic size.
 
template<typename DataType >
void registerStdVectorOfVars (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, std::vector< DataType > &data, const bool loadingMode=false)
 Method for registering a std::vector<DataType> of primitive DataType (int, double, ...)
 
template<typename DataType >
void registerStdVectorOfSerializablesOfConstSize (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, std::vector< DataType > &data, const bool loadingMode=false)
 Method for registering a std::vector<DataType> of constant-sized Serializable
 
template<typename DataType >
void registerStdVectorOfSerializables (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, std::vector< DataType > &data, const bool loadingMode=false)
 Method for registering a std::vector<DataType> of dynamic-sized DataType
 
template<typename DataTypeKey , typename DataTypeValue >
void registerMap (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, std::map< DataTypeKey, DataTypeValue > &data, const bool loadingMode=false)
 Method for registering a std::map<DataTypeKey, DataTypeValue> of fixed-sized types (i.e. int, double)
 
size_t addSizeToBuffer (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, size_t &sizeBufferIndex, bool *&dataPtr, const size_t data) const
 Add a size_t to the sizeBuffer in the n-th util::round and return that size_t in all successive rounds.
 
- Protected Member Functions inherited from olb::Serializable
template<typename DataType >
void registerVar (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, bool *&dataPtr, const DataType &data, const size_t arrayLength=1) const
 Register primitive data types (int, double, ...) or arrays of those.
 
template<typename DataType >
void registerSerializableOfConstSize (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, bool *&dataPtr, DataType &data, const bool loadingMode=false)
 Register Serializable object of constant size.
 
template<typename DataType >
void registerSerializablesOfConstSize (const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, bool *&dataPtr, DataType *data, const size_t arrayLength, const bool loadingMode=false)
 Register an array of Serializable objects of constant size.
 
- Protected Attributes inherited from olb::SuperStructure< T, DESCRIPTOR::d >
CuboidGeometry< T, D > & _cuboidGeometry
 The grid structure is stored here.
 
LoadBalancer< T > & _loadBalancer
 Distribution of the cuboids of the cuboid structure.
 
int _overlap
 Size of ghost cell layer (must be greater than 1 and greater_overlapBC, default =1)
 
OstreamManager clout
 class specific output stream
 
- Protected Attributes inherited from olb::BufferSerializable
std::vector< bool * > _dataBuffer
 Data buffer for data that has to be buffered between two getBlock() iterations.
 
std::vector< size_t > _sizeBuffer
 std::vector of integer buffers (e.g. for std::vector size) to be buffered for the whole iteration process
 

Detailed Description

template<typename T, typename DESCRIPTOR>
class olb::SuperLattice< T, DESCRIPTOR >

Super class maintaining block lattices for a cuboid decomposition.

Definition at line 45 of file superLattice.h.

Member Typedef Documentation

◆ block_t

template<typename T , typename DESCRIPTOR >
using olb::SuperLattice< T, DESCRIPTOR >::block_t = ConcretizableBlockLattice<T,DESCRIPTOR>

Definition at line 77 of file superLattice.h.

◆ descriptor_t

template<typename T , typename DESCRIPTOR >
using olb::SuperLattice< T, DESCRIPTOR >::descriptor_t = DESCRIPTOR

Descriptor / discrete velocity set of the lattice.

Definition at line 75 of file superLattice.h.

◆ value_t

template<typename T , typename DESCRIPTOR >
using olb::SuperLattice< T, DESCRIPTOR >::value_t = T

Base value type of the lattice.

Definition at line 73 of file superLattice.h.

Constructor & Destructor Documentation

◆ SuperLattice() [1/2]

template<typename T , typename DESCRIPTOR >
olb::SuperLattice< T, DESCRIPTOR >::SuperLattice ( SuperGeometry< T, DESCRIPTOR::d > & superGeometry)

Construct lattice for the cuboid decomposition of superGeometry.

Definition at line 101 of file superLattice.hh.

102 : SuperStructure<T,DESCRIPTOR::d>(superGeometry.getCuboidGeometry(),
103 superGeometry.getLoadBalancer(),
104 superGeometry.getOverlap()),
105 _statistics()
106{
107 using namespace stage;
108
109 auto& load = this->getLoadBalancer();
110
111 for (int iC = 0; iC < load.size(); ++iC) {
112 auto& cuboid = this->_cuboidGeometry.get(load.glob(iC));
113 #ifdef PLATFORM_GPU_CUDA
114 if (load.platform(iC) == Platform::GPU_CUDA) {
115 if (gpu::cuda::device::getCount() == 0) {
116 throw std::runtime_error("Load balancer requested GPU processing for cuboid "
117 + std::to_string(load.glob(iC))
118 + " but no CUDA device was found on rank "
119 + std::to_string(singleton::mpi().getRank()));
120
121 }
122 }
123 #endif
124 _block.emplace_back(constructUsingConcretePlatform<ConcretizableBlockLattice<T,DESCRIPTOR>>(
125 load.platform(iC), cuboid.getExtent(), this->getOverlap()));
126 }
127
128 {
129 auto& communicator = getCommunicator(PostCollide());
130 communicator.template requestField<descriptors::POPULATION>();
131 communicator.requestOverlap(1); // Required for inter-block propagation
132 communicator.exchangeRequests();
133 }
134
135 {
136 auto& communicator = getCommunicator(Full());
137 DESCRIPTOR::fields_t::for_each([&](auto field) {
138 communicator.template requestField<typename decltype(field)::type>();
139 });
140 // VTK output includes overlap of 1, some implicit dependencies for overlap of 2 exist
141 communicator.requestOverlap(std::min(2, this->getOverlap()));
142 communicator.exchangeRequests();
143 }
144
145 _statisticsEnabled = true;
146 _communicationNeeded = true;
147}
bool load(std::string fileName="", const bool enforceUint=false)
Load Serializable from file fileName
SuperCommunicator< T, SuperLattice > & getCommunicator(STAGE stage=STAGE())
Return communicator for given communication stage.
CuboidGeometry< T, D > & getCuboidGeometry()
Read and write access to cuboid geometry.
int getOverlap()
Read and write access to the overlap.
CuboidGeometry< T, D > & _cuboidGeometry
The grid structure is stored here.
LoadBalancer< T > & getLoadBalancer()
Read and write access to the load balancer.
int getRank() const
Returns the process ID.
int getCount()
Return number of available devices.
Definition device.hh:42
MpiManager & mpi()
@ GPU_CUDA
Vector CPU (AVX2 / AVX-512 collision)
CONCRETIZABLE::base_t * constructUsingConcretePlatform(Platform platform, ARGS &&... args)
Definition platform.h:133

References olb::SuperStructure< T, DESCRIPTOR::d >::_cuboidGeometry, olb::constructUsingConcretePlatform(), olb::SuperLattice< T, DESCRIPTOR >::getCommunicator(), olb::gpu::cuda::device::getCount(), olb::SuperStructure< T, DESCRIPTOR::d >::getLoadBalancer(), olb::SuperStructure< T, DESCRIPTOR::d >::getOverlap(), olb::singleton::MpiManager::getRank(), olb::GPU_CUDA, olb::Serializable::load(), and olb::singleton::mpi().

+ Here is the call graph for this function:

◆ SuperLattice() [2/2]

template<typename T , typename DESCRIPTOR >
olb::SuperLattice< T, DESCRIPTOR >::SuperLattice ( const SuperLattice< T, DESCRIPTOR > & )
delete

◆ ~SuperLattice()

template<typename T , typename DESCRIPTOR >
olb::SuperLattice< T, DESCRIPTOR >::~SuperLattice ( )
default

Member Function Documentation

◆ addCustomTask()

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addCustomTask ( std::function< void()> f)
inline

Schedules f for execution during every invocation of STAGE.

Definition at line 428 of file superLattice.h.

428 {
429 _customTasks[typeid(STAGE)].emplace_back(f);
430 }

◆ addLatticeCoupling() [1/6]

template<typename T , typename DESCRIPTOR >
template<typename... PARTNER_DESCRIPTORS>
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
PARTNER_DESCRIPTORS &... partnerLattices )

Adds a coupling generator for a multiple partner superLattice (legacy)

Definition at line 742 of file superLattice.hh.

745{
746 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
747 std::vector<BlockStructureD<DESCRIPTOR::d>*> partners{&partnerLattices.getBlock(iC)...};
748 getBlock(iC).addLatticeCoupling(indicator->getBlockIndicatorF(iC),
749 lcGen,
750 partners);
751 }
752}
BlockLattice< T, DESCRIPTOR > & getBlock(int locC)
Return BlockLattice with local index locC.
LoadBalancer< T > & _loadBalancer
Distribution of the cuboids of the cuboid structure.

◆ addLatticeCoupling() [2/6]

template<typename T , typename DESCRIPTOR >
template<typename PARTNER_DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices )

Adds a coupling generator for a vector of partner superLattice (legacy)

Definition at line 698 of file superLattice.hh.

701{
702 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
703 std::vector<BlockStructureD<DESCRIPTOR::d>*> partners;
704 for (auto partnerLattice: partnerLattices) {
705 partners.push_back(&partnerLattice->getBlock(iC));
706 }
707 getBlock(iC).addLatticeCoupling(indicator->getBlockIndicatorF(iC),
708 lcGen,
709 partners);
710 }
711}

◆ addLatticeCoupling() [3/6]

template<typename T , typename DESCRIPTOR >
template<typename... PARTNER_DESCRIPTORS>
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
PARTNER_DESCRIPTORS &... partnerLattices )

Adds a coupling generator for a multiple partner superLattice (legacy)

Definition at line 725 of file superLattice.hh.

728{
729 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
730 std::vector<BlockStructureD<DESCRIPTOR::d>*> partners{&partnerLattices.getBlock(iC)...};
731 int overlap = getBlock(iC).getPadding();
732 LatticeCouplingGenerator<T,DESCRIPTOR> *extractedLcGen = lcGen.clone();
733 extractedLcGen->reset(-overlap+1, getBlock(iC).getExtent()+overlap-2);
734 getBlock(iC).addLatticeCoupling(*extractedLcGen, partners);
735 delete extractedLcGen;
736 }
737 return;
738}

◆ addLatticeCoupling() [4/6]

template<typename T , typename DESCRIPTOR >
template<typename PARTNER_DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices )

Adds a coupling generator for a vector of partner superLattice (legacy)

Definition at line 678 of file superLattice.hh.

681{
682 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
683 std::vector<BlockStructureD<DESCRIPTOR::d>*> partners;
684 for (auto partnerLattice: partnerLattices) {
685 partners.push_back(&partnerLattice->getBlock(iC));
686 }
687 int overlap = getBlock(iC).getPadding();
688 LatticeCouplingGenerator<T,DESCRIPTOR> *extractedLcGen = lcGen.clone();
689 extractedLcGen->reset(-overlap+1, getBlock(iC).getExtent()+overlap-2);
690 getBlock(iC).addLatticeCoupling(*extractedLcGen, partners);
691 delete extractedLcGen;
692 }
693 return;
694}

◆ addLatticeCoupling() [5/6]

template<typename T , typename DESCRIPTOR >
template<typename... PARTNER_DESCRIPTORS>
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
PARTNER_DESCRIPTORS &... partnerLattices )

Adds a coupling generator for a multiple partner superLattice (legacy)

Definition at line 756 of file superLattice.hh.

759{
760 addLatticeCoupling(sGeometry.getMaterialIndicator(material),
761 lcGen, partnerLattices...);
762}
std::unique_ptr< SuperIndicatorF< T, D > > getMaterialIndicator(std::vector< int > &&materials)
Returns a material indicator using the given vector of materials.
void addLatticeCoupling(LatticeCouplingGenerator< T, DESCRIPTOR > const &lcGen, std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices)
Adds a coupling generator for a vector of partner superLattice (legacy)

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ addLatticeCoupling() [6/6]

template<typename T , typename DESCRIPTOR >
template<typename PARTNER_DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::addLatticeCoupling ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
LatticeCouplingGenerator< T, DESCRIPTOR > const & lcGen,
std::vector< SuperLattice< T, PARTNER_DESCRIPTOR > * > partnerLattices )

Adds a coupling generator for a vector of partner superLattice (legacy)

Definition at line 715 of file superLattice.hh.

718{
719 addLatticeCoupling(sGeometry.getMaterialIndicator(material),
720 lcGen, partnerLattices);
721}

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ addPostProcessor() [1/5]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addPostProcessor ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
PostProcessorGenerator< T, DESCRIPTOR > const & ppGen )

Add a non-local post-processing step (legacy)

Definition at line 635 of file superLattice.hh.

637{
638 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
639 getBlock(iC).template addPostProcessor<STAGE>(indicator->getBlockIndicatorF(iC),
640 ppGen);
641 }
642 _communicationNeeded = true;
643}

◆ addPostProcessor() [2/5]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addPostProcessor ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
PostProcessorPromise< T, DESCRIPTOR > && promise )

Add a non-local post-processing step (legacy)

Definition at line 647 of file superLattice.hh.

649{
650 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
651 getBlock(iC).addPostProcessor(typeid(STAGE),
652 indicator->getBlockIndicatorF(iC),
653 std::forward<decltype(promise)>(promise));
654 }
655}

◆ addPostProcessor() [3/5]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addPostProcessor ( PostProcessorGenerator< T, DESCRIPTOR > const & ppGen)

Add a non-local post-processing step (legacy)

Definition at line 621 of file superLattice.hh.

623{
624 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
625 int overlap = getBlock(iC).getPadding();
626 PostProcessorGenerator<T,DESCRIPTOR> *extractedPpGen = ppGen.clone();
627 extractedPpGen->reset(-overlap+1, getBlock(iC).getExtent()+overlap-2);
628 getBlock(iC).template addPostProcessor<STAGE>(*extractedPpGen);
629 delete extractedPpGen;
630 }
631}

References olb::LatticeStatistics< T >::reset().

+ Here is the call graph for this function:

◆ addPostProcessor() [4/5]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addPostProcessor ( PostProcessorPromise< T, DESCRIPTOR > && promise)

Add a non-local post-processing step (legacy)

Definition at line 659 of file superLattice.hh.

660{
661 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
662 getBlock(iC).addPostProcessor(typeid(STAGE),
663 std::forward<decltype(promise)>(promise));
664 }
665}

◆ addPostProcessor() [5/5]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::addPostProcessor ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
PostProcessorGenerator< T, DESCRIPTOR > const & ppGen )

Add a non-local post-processing step (legacy)

Definition at line 669 of file superLattice.hh.

672{
673 addPostProcessor<STAGE>(sGeometry.getMaterialIndicator(material), ppGen);
674}

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ collideAndStream()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::collideAndStream ( )

Core implementation of a single iteration of the collide and stream loop.

  1. Pre-collision communication (optional)
  1. Collide interior of all local block lattices
  2. Post-collision communicate for inter-block propagation
  3. Block local streaming
  4. Post-stream communication for boundary conditions / post processors (optional)
  5. Execute default post processors on all local block lattices
  6. Execute manually scheduled post processing callables (optional)
  7. Post-post-process communication (optional)
  8. Reset lattice statistics and mark overlap state as unclean

Definition at line 512 of file superLattice.hh.

513{
514 using namespace stage;
515
516 waitForBackgroundTasks(PreCollide());
517 auto& load = this->_loadBalancer;
518
519 if (_statisticsEnabled) {
520 setParameter<statistics::AVERAGE_RHO>(getStatistics().getAverageRho());
521 }
522
523 // Optional pre processing stage
524 executePostProcessors(PreCollide());
525
526 #ifdef PARALLEL_MODE_OMP
527 #pragma omp taskloop
528 #endif
529 for (int iC = 0; iC < load.size(); ++iC) {
530 _block[iC]->collide();
531 }
532
533 // Communicate propagation overlap, optional post processing
534 executePostProcessors(PostCollide());
535
536 // Block-local propagation
537 for (int iC = 0; iC < load.size(); ++iC) {
538 _block[iC]->stream();
539 }
540
541 // Communicate (default) post processor neighborhood and apply them
542 executePostProcessors(PostStream());
543
544 // Execute custom tasks (arbitrary callables)
545 // (used for multi-stage models such as free surface)
546 executeCustomTasks(PostStream());
547
548 // Final communication stage (e.g. for external coupling)
549 getCommunicator(PostPostProcess()).communicate();
550
551 if (_statisticsEnabled) {
552 collectStatistics();
553 }
554 _communicationNeeded = true;
555}
void waitForBackgroundTasks(STAGE stage=STAGE{})
Block until all background tasks scheduled for STAGE are completed.
LatticeStatistics< T > & getStatistics()
Return a handle to the LatticeStatistics object.
void executeCustomTasks(STAGE stage=STAGE())
Executes custom tasks assigned to STAGE.
void executePostProcessors(STAGE stage=STAGE())
Executes post processors for STAGE.

◆ communicate()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::communicate ( )
overridevirtual

Perform full overlap communication if needed.

Reimplemented from olb::SuperStructure< T, DESCRIPTOR::d >.

Definition at line 590 of file superLattice.hh.

591{
592 if (_communicationNeeded) {
593 getCommunicator(stage::Full()).communicate();
594 _communicationNeeded = false;
595 }
596}
+ Here is the caller graph for this function:

◆ defineDynamics() [1/6]

template<typename T , typename DESCRIPTOR >
template<template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator)

Set dynamics of indicated cells to DYNAMICS.

Definition at line 191 of file superLattice.hh.

192{
193 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
194 _block[iC]->template defineDynamics<DYNAMICS>(indicator->getBlockIndicatorF(iC));
195 }
196}

◆ defineDynamics() [2/6]

template<typename T , typename DESCRIPTOR >
template<template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator)

Set dynamics of indicated cells to DYNAMICS<T,DESCRIPTOR>

◆ defineDynamics() [3/6]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
Dynamics< T, DESCRIPTOR > * dynamics )

Define the dynamics on a domain described by an indicator (legacy)

Definition at line 222 of file superLattice.hh.

224{
225 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
226 _block[iC]->defineDynamics(indicator->getBlockIndicatorF(iC), dynamics);
227 }
228}

◆ defineDynamics() [4/6]

template<typename T , typename DESCRIPTOR >
template<template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material )

Set dynamics of material cells to DYNAMICS.

Definition at line 200 of file superLattice.hh.

201{
202 defineDynamics<DYNAMICS>(sGeometry.getMaterialIndicator(material));
203}

◆ defineDynamics() [5/6]

template<typename T , typename DESCRIPTOR >
template<template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material )

Set dynamics of indicated cells to DYNAMICS<T,DESCRIPTOR>

◆ defineDynamics() [6/6]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineDynamics ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
Dynamics< T, DESCRIPTOR > * dynamics )

Define the dynamics on a domain with a particular material number (legacy)

Definition at line 231 of file superLattice.hh.

233{
234 defineDynamics(sGeometry.getMaterialIndicator(material), dynamics);
235}
void defineDynamics(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator)
Set dynamics of indicated cells to DYNAMICS.

◆ defineField() [1/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & field )

Defines a field on a domain described by an indicator.

Definition at line 366 of file superLattice.hh.

368{
369 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
370 _block[iC]->template defineField<FIELD>(indicator->getBlockIndicatorF(iC),
371 field);
372 }
373 _communicationNeeded = true;
374}

◆ defineField() [2/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
FunctorPtr< SuperF< DESCRIPTOR::d, T, T > > && field )

Define an external field on a domain described by an indicator.

Definition at line 336 of file superLattice.hh.

338{
339 if (field->getBlockFSize() == this->_loadBalancer.size()) {
340 for (int iC=0; iC < this->_loadBalancer.size(); ++iC) {
341 _block[iC]->template defineField<FIELD>(indicator->getBlockIndicatorF(iC),
342 field->getBlockF(iC));
343 }
344 }
345 else {
346 FieldD<T,DESCRIPTOR,FIELD> fieldTmp;
347 int coords[DESCRIPTOR::d+1];
348 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
349 _block[iC]->forSpatialLocations([&](LatticeR<DESCRIPTOR::d> loc) {
350 coords[0] = iC;
351 for (int j = 0; j < DESCRIPTOR::d; ++j) {
352 coords[j+1] = loc[j];
353 }
354 if (indicator(coords)) {
355 field(fieldTmp.data(), coords);
356 _block[iC]->get(loc).template setField<FIELD>(fieldTmp);
357 }
358 });
359 }
360 }
361 _communicationNeeded = true;
362}
std::enable_if_t< DESCRIPTOR::d==2, std::shared_ptr< SuperF2D< T > > > field(SuperLattice< T, DESCRIPTOR > &sLattice)
Returns external field functor.

◆ defineField() [3/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
IndicatorF2D< T > & indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & field )

Defines a field on a indicated domain.

Definition at line 395 of file superLattice.hh.

397{
398 SuperIndicatorFfromIndicatorF2D<T> indicatorF(indicator, sGeometry);
399 defineField<FIELD>(indicatorF, field);
400}

◆ defineField() [4/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
IndicatorF3D< T > & indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & field )

Definition at line 404 of file superLattice.hh.

406{
407 SuperIndicatorFfromIndicatorF3D<T> indicatorF(indicator, sGeometry);
408 defineField<FIELD>(indicatorF, field);
409}

◆ defineField() [5/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & field )

Defines a field on a domain with a particular material number.

Definition at line 386 of file superLattice.hh.

389{
390 defineField<FIELD>(sGeometry.getMaterialIndicator(material), field);
391}

◆ defineField() [6/6]

template<typename T , typename DESCRIPTOR >
template<typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::defineField ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
FunctorPtr< SuperF< DESCRIPTOR::d, T, T > > && field )

Define an external field on a domain with a particular material number.

Definition at line 378 of file superLattice.hh.

380{
381 defineField<FIELD>(sGeometry.getMaterialIndicator(material), std::forward<decltype(field)>(field));
382}

◆ definePopulations() [1/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::definePopulations ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & Pop )

Define a population on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
PopAnalytical functor, target dimension DESCRIPTOR::q

Definition at line 300 of file superLattice.hh.

302{
303 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
304 _block[iC]->definePopulations(indicator->getBlockIndicatorF(iC), Pop);
305 }
306 _communicationNeeded = true;
307}

◆ definePopulations() [2/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::definePopulations ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
SuperF< DESCRIPTOR::d, T, T > & Pop )
Parameters
indicatorIndicator describing the target domain
PopSuper functor, target dimension DESCRIPTOR::q

Definition at line 317 of file superLattice.hh.

319{
320 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
321 _block[iC]->definePopulations(indicator->getBlockIndicatorF(iC), Pop.getBlockF(iC));
322 }
323 _communicationNeeded = true;
324}

◆ definePopulations() [3/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::definePopulations ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & Pop )

Define a population on a domain with a particular material number.

Definition at line 310 of file superLattice.hh.

312{
313 definePopulations(sGeometry.getMaterialIndicator(material), Pop);
314}
void definePopulations(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &Pop)
Define a population on a domain described by an indicator.

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ definePopulations() [4/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::definePopulations ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
SuperF< DESCRIPTOR::d, T, T > & Pop )

Define a population on a domain with a particular material number.

Definition at line 327 of file superLattice.hh.

329{
330 definePopulations(sGeometry.getMaterialIndicator(material), Pop);
331}

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ defineRho() [1/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineRho ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & rho )

Define rho on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
rhoAnalytical functor

Definition at line 238 of file superLattice.hh.

240{
241 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
242 _block[iC]->defineRho(indicator->getBlockIndicatorF(iC),
243 rho);
244 }
245 _communicationNeeded = true;
246}

◆ defineRho() [2/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineRho ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & rho )

Define rho on a domain with a particular material number.

Definition at line 249 of file superLattice.hh.

251{
252 defineRho(sGeometry.getMaterialIndicator(material), rho);
253}
void defineRho(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&, AnalyticalF< DESCRIPTOR::d, T, T > &rho)
Define rho on a domain described by an indicator.

◆ defineRhoU() [1/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineRhoU ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Define rho and u on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
rhoAnalytical functor
uAnalytical functor

Definition at line 277 of file superLattice.hh.

280{
281 #ifdef PARALLEL_MODE_OMP
282 #pragma omp parallel for schedule(dynamic,1)
283 #endif
284 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
285 _block[iC]->defineRhoU(indicator->getBlockIndicatorF(iC),
286 rho, u);
287 }
288 _communicationNeeded = true;
289}

◆ defineRhoU() [2/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineRhoU ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Define rho and u on a domain with a particular material number.

Definition at line 292 of file superLattice.hh.

295{
296 defineRhoU(sGeometry.getMaterialIndicator(material), rho, u);
297}
void defineRhoU(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
Define rho and u on a domain described by an indicator.

◆ defineU() [1/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineU ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Define u on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
uAnalytical functor

Definition at line 256 of file superLattice.hh.

258{
259 #ifdef PARALLEL_MODE_OMP
260 #pragma omp parallel for schedule(dynamic,1)
261 #endif
262 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
263 _block[iC]->defineU(indicator->getBlockIndicatorF(iC),
264 u);
265 }
266 _communicationNeeded = true;
267}

◆ defineU() [2/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::defineU ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Define u on a domain with a particular material number.

Definition at line 270 of file superLattice.hh.

272{
273 defineU(sGeometry.getMaterialIndicator(material), u);
274}
void defineU(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &u)
Define u on a domain described by an indicator.

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ executeCoupling()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::executeCoupling ( )

Executes coupling with partner lattices (legacy)

  1. Pre-coupling communication (optional)
  2. Execute coupling post processors on all local block lattices
  3. Post-coupling communication (optional)

Definition at line 765 of file superLattice.hh.

766{
767 getCommunicator(stage::PreCoupling()).communicate();
768 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
769 _block[iC]->executeCoupling();
770 }
771 executeCustomTasks(stage::Coupling());
772 getCommunicator(stage::PostCoupling()).communicate();
773 _communicationNeeded = true;
774}

◆ executeCustomTasks()

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::executeCustomTasks ( STAGE stage = STAGE())

Executes custom tasks assigned to STAGE.

Definition at line 778 of file superLattice.hh.

779{
780 for (auto& f : _customTasks[typeid(STAGE)]) {
781 f();
782 }
783}

◆ executePostProcessors()

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::executePostProcessors ( STAGE stage = STAGE())

Executes post processors for STAGE.

  1. Pre-communication for STAGE
  2. Execute post processors for STAGE on all block lattices

Definition at line 559 of file superLattice.hh.

560{
561 #ifdef PLATFORM_GPU_CUDA
563 #endif
564
565 getCommunicator(stage).communicate();
566
567 auto& load = this->_loadBalancer;
568
569 #ifdef PARALLEL_MODE_OMP
570 #pragma omp taskloop
571 #endif
572 for (int iC = 0; iC < load.size(); ++iC) {
573 _block[iC]->template postProcess<STAGE>();
574 }
575}
void synchronize()
Synchronize device.
Definition device.hh:64

References olb::gpu::cuda::device::synchronize().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ forBlocksOnPlatform()

template<typename T , typename DESCRIPTOR >
template<Platform PLATFORM, typename F >
void olb::SuperLattice< T, DESCRIPTOR >::forBlocksOnPlatform ( F f)
inline

Apply f to every ConcreteBlockLattice of PLATFORM.

Definition at line 112 of file superLattice.h.

112 {
113 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
114 if (_block[iC]->getPlatform() == PLATFORM) {
115 f(static_cast<ConcreteBlockLattice<T,DESCRIPTOR,PLATFORM>&>(*_block[iC]));
116 }
117 }
118 };

References olb::SuperStructure< T, DESCRIPTOR::d >::_loadBalancer.

◆ get() [1/2]

template<typename T , typename DESCRIPTOR >
Cell< T, DESCRIPTOR > olb::SuperLattice< T, DESCRIPTOR >::get ( LatticeR< DESCRIPTOR::d+1 > latticeR)

Get local cell interface.

Definition at line 150 of file superLattice.hh.

151{
152#ifdef PARALLEL_MODE_MPI
153 if (this->_loadBalancer.isLocal(latticeR[0])) {
154 if constexpr (DESCRIPTOR::d == 3) {
155 return _block[this->_loadBalancer.loc(latticeR[0])]->get(latticeR[1], latticeR[2], latticeR[3]);
156 } else {
157 return _block[this->_loadBalancer.loc(latticeR[0])]->get(latticeR[1], latticeR[2]);
158 }
159 }
160 else {
161 throw std::domain_error("Cuboid iC must be locally available");
162 }
163#else
164 if constexpr (DESCRIPTOR::d == 3) {
165 return _block[this->_loadBalancer.loc(latticeR[0])]->get(latticeR[1], latticeR[2], latticeR[3]);
166 } else {
167 return _block[this->_loadBalancer.loc(latticeR[0])]->get(latticeR[1], latticeR[2]);
168 }
169#endif
170}

◆ get() [2/2]

template<typename T , typename DESCRIPTOR >
template<typename... R>
std::enable_if_t< sizeof...(R)==DESCRIPTOR::d+1, Cell< T, DESCRIPTOR > > olb::SuperLattice< T, DESCRIPTOR >::get ( R... latticeR)

Get local cell interface.

Definition at line 175 of file superLattice.hh.

176{
177 return get({latticeR...});
178}
Cell< T, DESCRIPTOR > get(LatticeR< DESCRIPTOR::d+1 > latticeR)
Get local cell interface.

◆ getBlock() [1/5]

template<typename T , typename DESCRIPTOR >
BlockLattice< T, DESCRIPTOR > & olb::SuperLattice< T, DESCRIPTOR >::getBlock ( int locC)
inline

Return BlockLattice with local index locC.

Definition at line 85 of file superLattice.h.

86 {
87 return *_block[locC];
88 };
+ Here is the caller graph for this function:

◆ getBlock() [2/5]

template<typename T , typename DESCRIPTOR >
template<typename BLOCK >
BLOCK & olb::SuperLattice< T, DESCRIPTOR >::getBlock ( int locC)
inline

Return locC-th block lattice casted to BLOCK.

Only safe to use under external knownledge that the locC-th block is of type BLOCK

Definition at line 94 of file superLattice.h.

95 {
96 return *dynamic_cast<BLOCK*>(_block[locC].get());
97 };

◆ getBlock() [3/5]

template<typename T , typename DESCRIPTOR >
const BlockLattice< T, DESCRIPTOR > & olb::SuperLattice< T, DESCRIPTOR >::getBlock ( int locC) const
inline

Return read-only BlockLattice with local index locC.

Definition at line 99 of file superLattice.h.

100 {
101 return *_block[locC];
102 };

◆ getBlock() [4/5]

template<typename T , typename DESCRIPTOR >
template<typename BLOCK >
const BLOCK & olb::SuperLattice< T, DESCRIPTOR >::getBlock ( int locIC) const
inline

Return locC-th block lattice casted to BLOCK (read-only)

Definition at line 105 of file superLattice.h.

106 {
107 return *dynamic_cast<const BLOCK*>(_block[locIC].get());
108 };

◆ getBlock() [5/5]

template<typename T , typename DESCRIPTOR >
bool * olb::SuperLattice< T, DESCRIPTOR >::getBlock ( std::size_t iBlock,
std::size_t & sizeBlock,
bool loadingMode )
overridevirtual

Return a pointer to the memory of the current block and its size for the serializable interface.

Implements olb::Serializable.

Definition at line 851 of file superLattice.hh.

852{
853 std::size_t currentBlock = 0;
854 bool* dataPtr = nullptr;
855
856 for (std::size_t iC=0; iC < _block.size(); ++iC) {
857 registerSerializableOfConstSize(iBlock, sizeBlock, currentBlock, dataPtr, getBlock(iC), loadingMode);
858 }
859
860 return dataPtr;
861}
void registerSerializableOfConstSize(const std::size_t iBlock, std::size_t &sizeBlock, std::size_t &currentBlock, bool *&dataPtr, DataType &data, const bool loadingMode=false)
Register Serializable object of constant size.
Definition serializer.h:239

◆ getCommunicator() [1/2]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
SuperCommunicator< T, SuperLattice< T, DESCRIPTOR > > & olb::SuperLattice< T, DESCRIPTOR >::getCommunicator ( STAGE stage)

Definition at line 579 of file superLattice.hh.

580{
581 auto iter = _communicator.find(typeid(STAGE));
582 if (iter == _communicator.end()) {
583 iter = std::get<0>(_communicator.emplace(typeid(STAGE),
584 std::make_unique<SuperCommunicator<T,SuperLattice>>(*this)));
585 }
586 return *std::get<1>(*iter);
587}

◆ getCommunicator() [2/2]

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
SuperCommunicator< T, SuperLattice > & olb::SuperLattice< T, DESCRIPTOR >::getCommunicator ( STAGE stage = STAGE())

Return communicator for given communication stage.

+ Here is the caller graph for this function:

◆ getNblock()

template<typename T , typename DESCRIPTOR >
std::size_t olb::SuperLattice< T, DESCRIPTOR >::getNblock ( ) const
overridevirtual

Number of data blocks for the serializable interface.

Implements olb::Serializable.

Definition at line 834 of file superLattice.hh.

835{
836 return std::accumulate(_block.begin(), _block.end(), size_t(0), [](std::size_t sum, auto& b) -> std::size_t {
837 return sum + b->getNblock();
838 });
839}

◆ getSerializableSize()

template<typename T , typename DESCRIPTOR >
std::size_t olb::SuperLattice< T, DESCRIPTOR >::getSerializableSize ( ) const
overridevirtual

Binary size for the serializer.

Implements olb::Serializable.

Definition at line 843 of file superLattice.hh.

844{
845 return std::accumulate(_block.begin(), _block.end(), size_t(0), [](std::size_t sum, auto& b) -> std::size_t {
846 return sum + b->getSerializableSize();
847 });
848}

◆ getStatistics() [1/2]

template<typename T , typename DESCRIPTOR >
LatticeStatistics< T > & olb::SuperLattice< T, DESCRIPTOR >::getStatistics ( )

Return a handle to the LatticeStatistics object.

Definition at line 608 of file superLattice.hh.

609{
610 return _statistics;
611}

◆ getStatistics() [2/2]

template<typename T , typename DESCRIPTOR >
LatticeStatistics< T > const & olb::SuperLattice< T, DESCRIPTOR >::getStatistics ( ) const

Return a constant handle to the LatticeStatistics object.

Definition at line 614 of file superLattice.hh.

615{
616 return _statistics;
617}

◆ iniEquilibrium() [1/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniEquilibrium ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Initialize by equilibrium on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
rhoAnalytical functor (global)
uAnalytical functor (global)

Definition at line 452 of file superLattice.hh.

455{
456 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
457 _block[iC]->iniEquilibrium(indicator->getBlockIndicatorF(iC), rho, u);
458 }
459 _communicationNeeded = true;
460}

◆ iniEquilibrium() [2/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniEquilibrium ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
SuperF< DESCRIPTOR::d, T, T > & u )

Definition at line 471 of file superLattice.hh.

474{
475 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
476 _block[iC]->iniEquilibrium(indicator->getBlockIndicatorF(iC), rho, u.getBlockF(iC));
477 }
478 _communicationNeeded = true;
479}

◆ iniEquilibrium() [3/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniEquilibrium ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u )

Initialize by equilibrium on a domain with a particular material number.

Definition at line 463 of file superLattice.hh.

466{
467 iniEquilibrium(sGeometry.getMaterialIndicator(material), rho, u);
468}
void iniEquilibrium(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u)
Initialize by equilibrium on a domain described by an indicator.

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ iniEquilibrium() [4/4]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniEquilibrium ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
SuperF< DESCRIPTOR::d, T, T > & u )

Definition at line 482 of file superLattice.hh.

485{
486 iniEquilibrium(sGeometry.getMaterialIndicator(material), rho, u);
487}

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ iniRegularized() [1/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniRegularized ( FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > && indicator,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u,
AnalyticalF< DESCRIPTOR::d, T, T > & pi )

Initialize by non- and equilibrium on a domain described by an indicator.

Parameters
indicatorIndicator describing the target domain
rhoAnalytical functor (global)
uAnalytical functor (global)
piAnalytical functor (global)

Definition at line 491 of file superLattice.hh.

495{
496 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
497 _block[iC]->iniRegularized(indicator->getBlockIndicatorF(iC), rho, u, pi);
498 }
499 _communicationNeeded = true;
500}

◆ iniRegularized() [2/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::iniRegularized ( SuperGeometry< T, DESCRIPTOR::d > & sGeometry,
int material,
AnalyticalF< DESCRIPTOR::d, T, T > & rho,
AnalyticalF< DESCRIPTOR::d, T, T > & u,
AnalyticalF< DESCRIPTOR::d, T, T > & pi )

Initialize by non- and equilibrium on a domain with a particular material number.

Definition at line 503 of file superLattice.hh.

507{
508 iniRegularized(sGeometry.getMaterialIndicator(material), rho, u, pi);
509}
void iniRegularized(FunctorPtr< SuperIndicatorF< T, DESCRIPTOR::d > > &&indicator, AnalyticalF< DESCRIPTOR::d, T, T > &rho, AnalyticalF< DESCRIPTOR::d, T, T > &u, AnalyticalF< DESCRIPTOR::d, T, T > &pi)
Initialize by non- and equilibrium on a domain described by an indicator.

References olb::SuperGeometry< T, D >::getMaterialIndicator().

+ Here is the call graph for this function:

◆ initialize()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::initialize ( )

Initialize lattice to be ready for simulation.

Calls default post processors (PostStream) once and switches processing context to simulation.

Definition at line 181 of file superLattice.hh.

182{
183 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
184 _block[iC]->initialize();
185 }
186 _communicationNeeded = true;
187}

◆ postLoad()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::postLoad ( )
overridevirtual

Reimplemented from olb::Serializable.

Definition at line 864 of file superLattice.hh.

865{
866 for (int iC=0; iC < this->_loadBalancer.size(); ++iC) {
867 _block[iC]->postLoad();
868 }
869}

◆ scheduleBackgroundOutput()

template<typename T , typename DESCRIPTOR >
template<typename F >
void olb::SuperLattice< T, DESCRIPTOR >::scheduleBackgroundOutput ( F && f)

Schedule F(iC) for one-off background output of block data.

Definition at line 794 of file superLattice.hh.

795{
796 auto& load = this->getLoadBalancer();
797 for (int iC=0; iC < load.size(); ++iC) {
798 if (isPlatformCPU(load.platform(iC))) {
799 scheduleBackgroundTask<stage::PreCollide>(std::bind(f, iC));
800 }
801 }
802 for (int iC=0; iC < load.size(); ++iC) {
803 if (load.platform(iC) == Platform::GPU_CUDA) {
804 scheduleBackgroundTask<stage::PreContextSwitchTo<ProcessingContext::Evaluation>>(std::bind(f, iC));
805 }
806 }
807}
constexpr bool isPlatformCPU(Platform platform)
Returns true if platform is equal to Platform::CPU_*.
Definition platform.h:154

References olb::GPU_CUDA, and olb::isPlatformCPU().

+ Here is the call graph for this function:

◆ scheduleBackgroundOutputVTK()

template<typename T , typename DESCRIPTOR >
template<typename CONTEXT >
void olb::SuperLattice< T, DESCRIPTOR >::scheduleBackgroundOutputVTK ( CONTEXT && vtkContext)

Schedule one-off background output of given VTK CONTEXT.

CONTEXT must accept tasks as a callable task for which it provides SuperVTMwriter(2,3)D and iT arguments.

Definition at line 811 of file superLattice.hh.

812{
813 vtkContext([](auto& writer, std::size_t iT) {
814 writer.writePVD(iT);
815 });
816 scheduleBackgroundOutput([vtkContext](int iC) {
817 vtkContext([&](auto& writer, std::size_t iT) {
818 writer.writeVTI(iT, iC);
819 });
820 });
821}
void scheduleBackgroundOutput(F &&f)
Schedule F(iC) for one-off background output of block data.

◆ scheduleBackgroundTask()

template<typename T , typename DESCRIPTOR >
template<typename STAGE , typename F >
void olb::SuperLattice< T, DESCRIPTOR >::scheduleBackgroundTask ( F && f)

Schedule F for one-off background execution to be finished by STAGE.

Definition at line 787 of file superLattice.hh.

788{
789 _backgroundTasks[typeid(STAGE)].emplace_back(singleton::pool().schedule(f));
790}
std::future< R > schedule(F &&f)
Schedule F and return future of its return value.
Definition threadPool.h:145
ThreadPool & pool()
Definition olbInit.cpp:37

References olb::singleton::pool(), and olb::ThreadPool::schedule().

+ Here is the call graph for this function:

◆ setParameter() [1/2]

template<typename T , typename DESCRIPTOR >
template<typename PARAMETER , Platform PLATFORM, typename FIELD >
void olb::SuperLattice< T, DESCRIPTOR >::setParameter ( FieldArrayD< T, DESCRIPTOR, PLATFORM, FIELD > & fieldArray)

Set PARAMETER to column pointers of given field array.

PARAMETER must store size(FIELD) pointers to the value type of FIELD.

Useful for providing dynamically sized field arrays as parameters to certain post processors.

Definition at line 439 of file superLattice.hh.

440{
441 static_assert(DESCRIPTOR::template size<PARAMETER>() == DESCRIPTOR::template size<FIELD>(),
442 "PARAMETER size must equal FIELD size");
443 static_assert(std::is_same_v<typename PARAMETER::template value_type<T>,
444 typename FIELD::template value_type<T>*>,
445 "PARAMETER must store pointers to FIELD components");
446 for (int iC=0; iC < this->getLoadBalancer().size(); ++iC) {
447 _block[iC]->template setParameter<PARAMETER>(fieldArray);
448 }
449}

◆ setParameter() [2/2]

template<typename T , typename DESCRIPTOR >
template<typename PARAMETER >
void olb::SuperLattice< T, DESCRIPTOR >::setParameter ( FieldD< T, DESCRIPTOR, PARAMETER > field)

Update PARAMETER in all dynamics and post processors.

e.g. SuperLattice::setParameter<OMEGA>(0.6) sets the parameter OMEGA to 0.6 in all operators that declare it as one of their parameters.

Definition at line 413 of file superLattice.hh.

414{
415 for (int iC=0; iC < this->getLoadBalancer().size(); ++iC) {
416 _block[iC]->template setParameter<PARAMETER>(field);
417 }
418}

◆ setParameterOfDynamics() [1/2]

template<typename T , typename DESCRIPTOR >
template<typename PARAMETER , template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::setParameterOfDynamics ( FieldD< T, DESCRIPTOR, PARAMETER > && field)

Update PARAMETER in DYNAMICS.

Definition at line 422 of file superLattice.hh.

423{
424 for (int iC=0; iC < this->getLoadBalancer().size(); ++iC) {
425 _block[iC]->template getData<OperatorParameters<DYNAMICS>>().template set<PARAMETER>(
426 std::forward<decltype(field)>(field));
427 }
428}

◆ setParameterOfDynamics() [2/2]

template<typename T , typename DESCRIPTOR >
template<typename PARAMETER , template< typename... > typename DYNAMICS>
void olb::SuperLattice< T, DESCRIPTOR >::setParameterOfDynamics ( FieldD< T, DESCRIPTOR, PARAMETER > && field)

Update PARAMETER in DYNAMICS<T,DESCRIPTOR>

◆ setProcessingContext() [1/2]

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::setProcessingContext ( ProcessingContext context)
inline

Set processing context of block lattices.

Used to sync data between device and host

Definition at line 124 of file superLattice.h.

125 {
126 // Communicate overlap prior to evaluation
127 if (context == ProcessingContext::Evaluation) {
128 communicate();
129 waitForBackgroundTasks<stage::PreContextSwitchTo<ProcessingContext::Evaluation>>();
130 }
131 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
132 _block[iC]->setProcessingContext(context);
133 }
134 };
void communicate() override
Perform full overlap communication if needed.

References olb::SuperStructure< T, DESCRIPTOR::d >::_loadBalancer, olb::SuperLattice< T, DESCRIPTOR >::communicate(), and olb::Evaluation.

+ Here is the call graph for this function:

◆ setProcessingContext() [2/2]

template<typename T , typename DESCRIPTOR >
template<typename FIELD_TYPE >
void olb::SuperLattice< T, DESCRIPTOR >::setProcessingContext ( ProcessingContext context)
inline

Set processing context of FIELD_TYPE in block lattices.

Used to sync data between device and host

Definition at line 141 of file superLattice.h.

142 {
143 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
144 if (_block[iC]->template hasData<FIELD_TYPE>()) {
145 _block[iC]->template getData<FIELD_TYPE>().setProcessingContext(context);
146 }
147 }
148 };

References olb::SuperStructure< T, DESCRIPTOR::d >::_loadBalancer.

◆ statisticsOff()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::statisticsOff ( )
inline

Switch Statistics off (default on)

This speeds up the execution time

Definition at line 354 of file superLattice.h.

355 {
356 _statisticsEnabled = false;
357 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
358 _block[iC]->setStatisticsEnabled(false);
359 }
360 };

References olb::SuperStructure< T, DESCRIPTOR::d >::_loadBalancer.

◆ statisticsOn()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::statisticsOn ( )
inline

Switch Statistics on (default on)

Definition at line 343 of file superLattice.h.

344 {
345 _statisticsEnabled = true;
346 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
347 _block[iC]->setStatisticsEnabled(true);
348 }
349 };

References olb::SuperStructure< T, DESCRIPTOR::d >::_loadBalancer.

◆ stripeOffDensityOffset()

template<typename T , typename DESCRIPTOR >
void olb::SuperLattice< T, DESCRIPTOR >::stripeOffDensityOffset ( T offset)

Subtract constant offset from the density.

Definition at line 599 of file superLattice.hh.

600{
601 for (int iC = 0; iC < this->_loadBalancer.size(); ++iC) {
602 _block[iC]->stripeOffDensityOffset(offset);
603 }
604 _communicationNeeded = true;
605}

◆ waitForBackgroundTasks()

template<typename T , typename DESCRIPTOR >
template<typename STAGE >
void olb::SuperLattice< T, DESCRIPTOR >::waitForBackgroundTasks ( STAGE stage = STAGE{})

Block until all background tasks scheduled for STAGE are completed.

Definition at line 825 of file superLattice.hh.

826{
827 if (!_backgroundTasks[typeid(STAGE)].empty()) {
828 singleton::pool().waitFor(_backgroundTasks[typeid(STAGE)]);
829 _backgroundTasks[typeid(STAGE)].clear();
830 }
831}
void waitFor(std::vector< std::future< T > > &futures)
Blocks until all tasks producing the given futures are completed.
Definition threadPool.h:174

References olb::singleton::pool(), and olb::ThreadPool::waitFor().

+ Here is the call graph for this function:

Member Data Documentation

◆ d

template<typename T , typename DESCRIPTOR >
constexpr unsigned olb::SuperLattice< T, DESCRIPTOR >::d = DESCRIPTOR::d
staticconstexpr

Definition at line 70 of file superLattice.h.


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