Particle distribution for time and size discretization. More...

#include <ParticleDistribution.h>

Collaboration diagram for olb::ParticleDistribution< T, S >:

Public Member Functions
template<typename C >
	ParticleDistribution (FunctorPtr< AnalyticalF1D< T, S > > &&distribution, C numberOfParticles, C numberOfParticleSizeGroups, T minRadius, T maxRadius, T lengthConversion=1)

bool	calculateParticleArray (bool getExactParticleNumber=false)

template<typename DESCRIPTOR >
bool	calculateTimeArray (AnalyticalF1D< T, S > &timeDistribution, UnitConverter< T, DESCRIPTOR > &converter, T begin, T end)

void	printTimeArray ()
	prints out all time-steps a new particle should appear

template<typename C >
void	spawnSphericalParticles (SuperParticleSystem3D< T, Particle3D > &supParticleSystem, IndicatorF3D< S > &indicator, T density, C iT, bool loopAround=false, bool reshuffleAfterLoop=false)

template<typename C >
void	timeActivateParticles (SuperParticleSystem3D< T, Particle3D > &supParticleSystem, C iT, C idOfset=0)

void	deactivateParticles (SuperParticleSystem3D< T, Particle3D > &supParticleSystem, int beginID, int endID)
	Sets the particle state 'active' to 'false' for particles with ids from 'beginID' to 'endID'.

void	preSpawnSphericalParticles (SuperParticleSystem3D< T, Particle3D > &supParticleSystem, IndicatorF3D< S > &indicator, T density, int beginID=1, bool shuffle=false)

void	getParticleArray (std::vector< T > &particleArray)

T	nextParticleRadius (bool loopAround=false, bool reshuffleAfterLoop=false)

void	shuffleParticleArray ()
	shuffles particle array using random_shuffle

void	printSizeMatrix ()
	prints the calculation matrix to the console

void	printParticleArray ()
	prints the particle Radii to the console

	~ParticleDistribution ()
	destructor

Detailed Description

template<typename T, typename S>
class olb::ParticleDistribution< T, S >

Particle distribution for time and size discretization.

size discretization This class uses a given particle distribution ('distribution') and discretizes it in 'numberOfParticleSizeGroups' Radius-size categories, ranging from 'minRadius' to 'maxRadius'. The particles are stored in '_particleArray', where every entry equals the particle Radius (in meter) of one particle according to the distribution. To get the exact No of particles use 'calculateParticleArray(true)'. WARNING: leaving getExactParticleNumber = false may result in fewer or more particles than specified. BUT when putting getExactParticleNumber = true it may not converge. If the particles should appear random according to the distribution use 'shuffleParticleArray'. *

time discretization To calculate when a new particle shall be spawned, use 'calculateTimeArray' with the corresponding time distribution and the time interval ranging from 'begin' to 'end'. Every time-step ('iT') when a new particle should be spawned is stored in '_timeArray'. 'spawnSphericalParticles' can be called every time-step to automatically spawn a spherical particle according to the time and size discretization. Use 'loopAround = true', to start at the beginning of '_particleArray' after reaching its end and 'reshuffleAfterLoop = true' to reshuffle it at the end

time activation When dealing with large number of particles adding them will get more time consuming because of the needed communication between sub-grids of all threads. To counteract this, the particles can be added as passive particles in bulk with 'preSpawnSphericalParticles' and calling 'timeActivateParticles' will activate the particles according to the time distribution.

Definition at line 85 of file ParticleDistribution.h.

Constructor & Destructor Documentation

◆ ParticleDistribution()

template<typename T , typename S >

template<typename C >

olb::ParticleDistribution< T, S >::ParticleDistribution	(	FunctorPtr< AnalyticalF1D< T, S > > &&	distribution,
		C	numberOfParticles,
		C	numberOfParticleSizeGroups,
		T	minRadius,
		T	maxRadius,
		T	lengthConversion = 1 )

Definition at line 60 of file ParticleDistribution.hh.

                                                                                                                                                                                                  :
    _distribution(std::move(distribution)), _numberOfParticles(numberOfParticles), _numberOfParticleSizeGroups(numberOfParticleSizeGroups), _minRadius(minRadius), _maxRadius(maxRadius) {
    _sizeCalculation = new T * [_numberOfParticleSizeGroups];
    const T delta_d = (_maxRadius - _minRadius) / _numberOfParticleSizeGroups;
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      _sizeCalculation[i] = new T[4];
    }
    T output[1] = {};
    T input[1] = {};
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      // calculate size groups
      _sizeCalculation[i][0] = (_minRadius + delta_d * i) * lengthConversion;
      input[0] = _sizeCalculation[i][0]  / lengthConversion;
      // calculate concentration
      _distribution->operator()(output, input);
      _sizeCalculation[i][1] = output[0];
    }
    calculateParticleArray();
  };

References olb::ParticleDistribution< T, S >::calculateParticleArray().

Here is the call graph for this function:

◆ ~ParticleDistribution()

template<typename T , typename S >

olb::ParticleDistribution< T, S >::~ParticleDistribution ( )

destructor

Definition at line 331 of file ParticleDistribution.hh.

                                                    {
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      delete[] _sizeCalculation[i];
    }
    delete[] _sizeCalculation;
  };

Member Function Documentation

◆ calculateParticleArray()

template<typename T , typename S >

bool olb::ParticleDistribution< T, S >::calculateParticleArray ( bool getExactParticleNumber = false )

Definition at line 87 of file ParticleDistribution.hh.

                                                                                     {
    T sum_np = 0;
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      // calculate total particle conc.
      sum_np += _sizeCalculation[i][1];
    }
 
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      // calculate fraction
      _sizeCalculation[i][2] = _sizeCalculation[i][1] / sum_np;
      // calculate number of whole particles
      // ofset is added to correct difference in _numberOfParticles to _particleArray..size() from rounding errors
      _sizeCalculation[i][3] = round(_sizeCalculation[i][2] * (_numberOfParticles + _particleOffset));
    }
 
    long long p_count = 0;
    if (getExactParticleNumber) {
      for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
        p_count += _sizeCalculation[i][3];
      }
    }
    if (p_count == _numberOfParticles || getExactParticleNumber == false) {
      p_count = 0;
      for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
        if (_sizeCalculation[i][3] > 0) {
          for (long long j = 1; j <= _sizeCalculation[i][3]; j++) {
            // write particle radii to array, where every entry represents 1 particle
            _particleArray.push_back(_sizeCalculation[i][0]);
            p_count++;
          }
        }
      }
      return true;
    }
    else {
      _particleArray.clear();
      // calculate offset and calls function again with offset for convergence, choose
      // different (larger) divisor if function does not converge
      _particleOffset += (_numberOfParticles - p_count)/7.;
 
      return this->calculateParticleArray(true);
    }
  };

Here is the caller graph for this function:

◆ calculateTimeArray()

template<typename T , typename S >

template<typename DESCRIPTOR >

bool olb::ParticleDistribution< T, S >::calculateTimeArray	(	AnalyticalF1D< T, S > &	timeDistribution,
		UnitConverter< T, DESCRIPTOR > &	converter,
		T	begin,
		T	end )

Definition at line 200 of file ParticleDistribution.hh.

                                                                                                                                                    {
    _timeArray.clear();
    T input[1] = {  };
    T nSum = 0;
    T output[1] = { };
    T integral = 0;
    // calculate discrete integral
    for (std::size_t iT = converter.getLatticeTime(begin); iT <= converter.getLatticeTime(end); iT++) {
      input[0] = converter.getPhysTime(iT);
      timeDistribution(output, input);
 
      integral += output[0];
    }
    // evaluate timeDistribution for every time-step insde the time intervall
    for (std::size_t iT = converter.getLatticeTime(begin); iT <= converter.getLatticeTime(end); iT++) {
      input[0] = converter.getPhysTime(iT);
      timeDistribution(output, input);
      output[0] *= ((_numberOfParticles + _particleOffsetTime) / integral);
      if (output[0] > 0) {
        // adding up particle parts until a complete particle can be added
        nSum += output[0];
        if (nSum >= 1.) {
          for (int i = 1; i <= nSum; ++i) {
            // add the value of the time-step for one particle to _timeArray
            _timeArray.push_back(iT);
          }
          nSum = 0;
        }
      }
    }
    // check if all calculated particles match total number of particles
    if ((long long)(_timeArray.size()) == _numberOfParticles) {
      return true;
    }
    else {
      // calculate offset and calls function again with offset for convergence, choose
      // different (larger) divisor if function does not converge
      _particleOffsetTime += (_numberOfParticles - _timeArray.size())/7.;
      calculateTimeArray(timeDistribution, converter, begin, end);
      return false;
    }
  };

References olb::UnitConverter< T, DESCRIPTOR >::getLatticeTime(), and olb::UnitConverter< T, DESCRIPTOR >::getPhysTime().

Here is the call graph for this function:

◆ deactivateParticles()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::deactivateParticles	(	SuperParticleSystem3D< T, Particle3D > &	supParticleSystem,
		int	beginID,
		int	endID )

Sets the particle state 'active' to 'false' for particles with ids from 'beginID' to 'endID'.

Definition at line 299 of file ParticleDistribution.hh.

                                                                                                                                      {
    // Searches all sub-particle systems
      for (int pSystems = 0; pSystems < supParticleSystem.numOfPSystems(); pSystems++) {
        // searches all particles in sub-particle systems
        for (int particles = 0; particles < supParticleSystem.getParticleSystems()[pSystems]->size(); particles++) {
          // checks id of particle in sub-particle system
          int currentPiD = supParticleSystem.getParticleSystems()[pSystems]->getParticlesPointer()[particles]->getID();
          if (currentPiD <= endID && currentPiD >= beginID) {
            // deactivates particles
            supParticleSystem.getParticleSystems()[pSystems]->getParticlesPointer()[particles]->setActive(false);
          }
        }
      }
  };

References olb::SuperParticleSystem3D< T, PARTICLETYPE >::getParticleSystems(), and olb::SuperParticleSystem3D< T, PARTICLETYPE >::numOfPSystems().

Here is the call graph for this function:

◆ getParticleArray()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::getParticleArray ( std::vector< T > & particleArray )

Definition at line 131 of file ParticleDistribution.hh.

                                                                               {
    particleArray = _particleArray;
  };

◆ nextParticleRadius()

template<typename T , typename S >

T olb::ParticleDistribution< T, S >::nextParticleRadius	(	bool	loopAround = false,
		bool	reshuffleAfterLoop = false )

Definition at line 138 of file ParticleDistribution.hh.

                                                                                         {
    if (_currentParticle == 0 && shuffleParticles) {
      shuffleParticleArray();
      T Radius = _particleArray[_currentParticle];
      _currentParticle++;
      return Radius;
    }
    else if (_currentParticle == (_numberOfParticles)) {
      if (loopAround) {
        T Radius = _particleArray[_currentParticle];
        _currentParticle = 0;
        if (shuffleParticles) {
          shuffleParticleArray();
        }
        return Radius;
      }
      else {
        return 0;
      }
    }
    else {
      T Radius = _particleArray[_currentParticle];
      _currentParticle++;
      return Radius;
    }
  };

◆ preSpawnSphericalParticles()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::preSpawnSphericalParticles	(	SuperParticleSystem3D< T, Particle3D > &	supParticleSystem,
		IndicatorF3D< S > &	indicator,
		T	density,
		int	beginID = 1,
		bool	shuffle = false )

Definition at line 318 of file ParticleDistribution.hh.

                                                                                                                                                                                       {
    // call nextParticleRadius to get particle radius
    T p_radius = nextParticleRadius(false, shuffle);
    do {
      // adds a spherical particle with specified density and '_currentParticle' as its ID
      supParticleSystem.addTracerParticle(indicator, (T)(_currentParticle), 4. / 3. * M_PI * util::pow(p_radius, 3) * density, p_radius, 1);
      p_radius = nextParticleRadius(false, shuffle);
    } while (p_radius > 0);
    _currentParticle = 0;
    deactivateParticles(supParticleSystem, beginID, beginID + _numberOfParticles);
  };

References olb::SuperParticleSystem3D< T, PARTICLETYPE >::addTracerParticle(), M_PI, and olb::util::pow().

Here is the call graph for this function:

◆ printParticleArray()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::printParticleArray ( )

prints the particle Radii to the console

Definition at line 186 of file ParticleDistribution.hh.

                                                      {
    OstreamManager clout(std::cout, "printParticleArray");
    clout << std::setw(7) << "p_number" << std::setw(12) << "radius" << std::endl;
    for (long long i = 0; i < _particleArray.size(); i++) {
      clout << std::setw(7) << i << std::setw(12) << _particleArray[i] << std::endl;
    }
  };

◆ printSizeMatrix()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::printSizeMatrix ( )

prints the calculation matrix to the console

Definition at line 172 of file ParticleDistribution.hh.

                                                   {
    OstreamManager clout(std::cout, "printSizeMatrix");
    clout << std::setw(7) << "GroupNo" << std::setw(12) << "Radius" << std::setw(12) << "concentr." << std::setw(12) << "partiton" << std::setw(12) << "p_number" << std::endl;
    for (long long i = 0; i < _numberOfParticleSizeGroups; i++) {
      clout << std::setw(7) << i + 1;
      for (long long j = 0; j <= 3; j++) {
        clout << std::setw(12) << _sizeCalculation[i][j];
      }
      clout << std::endl;
    }
  };

◆ printTimeArray()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::printTimeArray ( )

prints out all time-steps a new particle should appear

Definition at line 245 of file ParticleDistribution.hh.

                                                  {
    OstreamManager clout(std::cout, "printTimeArray");
    clout << std::setw(7) << "p_number" << std::setw(12) << "time" << std::endl;
    for (long long i = 0; i < _timeArray.size(); i++) {
      clout << std::setw(7) << i << std::setw(12) << _timeArray[i] << std::endl;
    }
  };

◆ shuffleParticleArray()

template<typename T , typename S >

void olb::ParticleDistribution< T, S >::shuffleParticleArray ( )

shuffles particle array using random_shuffle

Definition at line 166 of file ParticleDistribution.hh.

                                                        {
    random_shuffle(&_particleArray[0], &_particleArray[_particleArray.size() - 1]);
  };

◆ spawnSphericalParticles()

template<typename T , typename S >

template<typename C >

void olb::ParticleDistribution< T, S >::spawnSphericalParticles	(	SuperParticleSystem3D< T, Particle3D > &	supParticleSystem,
		IndicatorF3D< S > &	indicator,
		T	density,
		C	iT,
		bool	loopAround = false,
		bool	reshuffleAfterLoop = false )

Definition at line 260 of file ParticleDistribution.hh.

                                                                                                                                                                                                       {
    while (_timeArray[_currentParticle] == iT) {
      T p_radius = nextParticleRadius(loopAround, shuffleParticles);
      if (p_radius > 0) {
        supParticleSystem.addParticle(indicator, 4. / 3. * M_PI * util::pow(p_radius, 3) * density, p_radius, 1);
      }
    }
  };

References olb::SuperParticleSystem3D< T, PARTICLETYPE >::addParticle(), M_PI, and olb::util::pow().

Here is the call graph for this function:

◆ timeActivateParticles()

template<typename T , typename S >

template<typename C >

void olb::ParticleDistribution< T, S >::timeActivateParticles	(	SuperParticleSystem3D< T, Particle3D > &	supParticleSystem,
		C	iT,
		C	idOfset = 0 )

Definition at line 274 of file ParticleDistribution.hh.

                                                                                                                                 {
    // check if a new particle has to be activated
    while (_timeArray[_currentParticle] == iT) {
      // Searches all sub-particle systems
      for (unsigned int pSystems = 0; pSystems < supParticleSystem.getParticleSystems().size(); pSystems++) {
        // checks id of particle in sub-particle system
        for (int particles = 0; particles < supParticleSystem.getParticleSystems()[pSystems]->size(); particles++) {
          unsigned int currentPiD = supParticleSystem.getParticleSystems()[pSystems]->getParticlesPointer()[particles]->getID();
          // checks id of particle in sub-particle system
          if (currentPiD == _currentParticle + 1 + idOfset) {
            // activates particles
            supParticleSystem.getParticleSystems()[pSystems]->getParticlesPointer()[particles]->setActive(true);
            }
        }
      }
      _currentParticle++;
      // starts at beginning when last particle is reached
      if(_currentParticle == _numberOfParticles){
        _currentParticle = 0;
      }
    }
  };

References olb::SuperParticleSystem3D< T, PARTICLETYPE >::getParticleSystems().

Here is the call graph for this function:

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

src/particles/subgrid3DLegacyFramework/ParticleDistribution.h
src/particles/subgrid3DLegacyFramework/ParticleDistribution.hh

Public Member Functions

Detailed Description

Constructor & Destructor Documentation

◆ ParticleDistribution()

◆ ~ParticleDistribution()

Member Function Documentation

◆ calculateParticleArray()

◆ calculateTimeArray()

◆ deactivateParticles()

◆ getParticleArray()

◆ nextParticleRadius()

◆ preSpawnSphericalParticles()

◆ printParticleArray()

◆ printSizeMatrix()

◆ printTimeArray()

◆ shuffleParticleArray()

◆ spawnSphericalParticles()

◆ timeActivateParticles()