revision: v0.14
PSD Class Reference

Contains a vector with radii and probabilities of a user defined particle size distribution (PSD) More...

#include <PSD.h>

Classes

class  RadiusAndProbability
 Class which stores radii and probabilities of a PSD. This class should be used as a vector<PSD::RadiusAndProbability>. More...
 

Public Types

enum  TYPE {
  TYPE::CUMULATIVE_NUMBER_DISTRIBUTION, TYPE::CUMULATIVE_LENGTH_DISTRIBUTION, TYPE::CUMULATIVE_VOLUME_DISTRIBUTION, TYPE::CUMULATIVE_AREA_DISTRIBUTION,
  TYPE::PROBABILITYDENSITY_NUMBER_DISTRIBUTION, TYPE::PROBABILITYDENSITY_LENGTH_DISTRIBUTION, TYPE::PROBABILITYDENSITY_AREA_DISTRIBUTION, TYPE::PROBABILITYDENSITY_VOLUME_DISTRIBUTION
}
 Enum class which stores the possible types of CDFs and PDFs. Particle size distributions can be represented by different probabilities based on number of particles, length of particles, surface area of particles or volume of particles. More...
 

Public Member Functions

 PSD ()
 Constructor; sets everything to 0 or default. More...
 
 PSD (const PSD &other)
 Copy constructor with deep copy. More...
 
 ~PSD ()
 Destructor; default destructor. More...
 
PSDcopy () const
 Creates a copy on the heap and returns a pointer. More...
 
void printPSD () const
 Prints radii and probabilities of the PSD vector. More...
 
Mdouble drawSample ()
 Draw a sample radius from a CUMULATIVE_NUMBER_DISTRIBUTION. More...
 
Mdouble insertManuallyByVolume (Mdouble volume)
 Draw sample radius manually per size class and check the volumeAllowed of each size class to insert the PSD as accurate as possible. More...
 
void validateCumulativeDistribution ()
 Validates if a CDF starts with zero and adds up to unity. More...
 
void validateProbabilityDensityDistribution ()
 Validates if the integral of the PDF equals to unity. More...
 
MERCURY_DEPRECATED void setPSDFromVector (std::vector< RadiusAndProbability > psd, TYPE PSDType)
 Deprecated version of reading in PSDs from a vector. More...
 
void setPSDFromCSV (const std::string &fileName, TYPE PSDType, bool headings=false, Mdouble unitScalingFactorRadii=1.0)
 read in the PSD vector with probabilities and radii saved in a .csv file. More...
 
void setDistributionUniform (Mdouble radMin, Mdouble radMax, int numberOfBins)
 create a PSD vector for a uniform distribution. More...
 
void setDistributionNormal (Mdouble mean, Mdouble standardDeviation, int numberOfBins)
 create a PSD vector for a normal distribution. More...
 
void convertProbabilityDensityToCumulative ()
 Converts a PDF to a CDF by integration. More...
 
void convertCumulativeToProbabilityDensity ()
 Converts a CDF to a PDF by derivation. More...
 
void convertProbabilityDensityToProbabilityDensityNumberDistribution (TYPE PDFType)
 convert any PDF to a PROBABILITYDENSITY_NUMBER_DISTRIBUTION. More...
 
void convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution ()
 convert a PROBABILITYDENSITY_NUMBER_DISTRIBUTION to a PROBABILITYDENSITY_VOLUME_DISTRIBUTION. More...
 
void convertCumulativeToCumulativeNumberDistribution (TYPE CDFType)
 convert any other CDF to a CUMULATIVE_NUMBER_DISTRIBUTION. More...
 
void cutoffCumulativeNumber (Mdouble quantileMin, Mdouble quantileMax, Mdouble minPolydispersity=0.1)
 cutoff the PSD at given quantiles. More...
 
void cutoffAndSqueezeCumulative (Mdouble quantileMin, Mdouble quantileMax, Mdouble squeeze, Mdouble minPolydispersity=0.1)
 cutoff the PSD at given quantiles and make it less polydisperse by squeezing it. More...
 
Mdouble getMinRadius () const
 Get smallest radius of the PSD. More...
 
Mdouble getMaxRadius () const
 Get largest radius of the PSD. More...
 
std::vector< RadiusAndProbabilitygetParticleSizeDistribution () const
 Get the PSD vector. More...
 
void setParticleSizeDistribution (std::vector< RadiusAndProbability >)
 set the PSD by a suitable vector. More...
 
int getInsertedParticleNumber () const
 Get the number of particles already inserted into the simulation. More...
 
void decrementNParticlesPerClass ()
 Decrement nParticlesPerClass_ counter. More...
 
void decrementVolumePerClass (Mdouble volume)
 Decrement volumePerClass_ counter. More...
 
Mdouble getNumberDx (Mdouble x) const
 Calculate a certain diameter (e.g. D10, D50, D90, etc.) from a percentile x of the number based PSD. More...
 
Mdouble getVolumeDx (Mdouble x) const
 Calculate a certain diameter (e.g. D10, D50, D90, etc.) from a percentile x of the volume based PSD. More...
 
Mdouble getRadiusByQuantile (Mdouble quantile) const
 Calculate the quantile of the PSD. More...
 
Mdouble getVolumetricMeanRadius () const
 get a volumetric mean radius of the PSD. More...
 
Mdouble getSizeRatio () const
 get the size ratio (width) of the PSD. More...
 
void cutHighSizeRatio ()
 Check if the size ratio is too high and cut it. More...
 
void computeRawMomenta ()
 compute raw momenta of the user defined PSD. More...
 
void computeCentralMomenta ()
 compute central momenta of the user defined PSD. More...
 
void computeStandardisedMomenta ()
 compute standardised momenta of the user defined PSD. More...
 
std::array< Mdouble, 6 > getMomenta () const
 get momenta of the user defined PSD. More...
 

Private Attributes

std::vector< RadiusAndProbabilityparticleSizeDistribution_
 
std::array< Mdouble, 6 > momenta_ {}
 
std::vector< intnParticlesPerClass_
 
std::vector< MdoublevolumePerClass_
 
int index_
 

Friends

bool operator< (const PSD::RadiusAndProbability &l, const PSD::RadiusAndProbability &r)
 determines if a certain value of the PSD vector is lower than another one. Used for std::lower_bound() More...
 
bool operator< (const PSD::RadiusAndProbability &l, Mdouble r)
 determines if a certain value of the PSD vector is lower than a double. More...
 
std::ostream & operator<< (std::ostream &os, PSD::RadiusAndProbability &p)
 Writes to output stream. More...
 
std::istream & operator>> (std::istream &is, PSD::RadiusAndProbability &p)
 Reads from input stream. More...
 
Mdouble operator== (PSD::RadiusAndProbability l, Mdouble r)
 Determines if a certain value of the PSD vector is equal to a double. More...
 

Detailed Description

Contains a vector with radii and probabilities of a user defined particle size distribution (PSD)

Stores radii and probabilities of a particle size distribution (PSD) in a vector of type PSD::RadiusAndProbability and converts them to a PSD which can be used by insertionBoundaries which insert particles into a simulation:

Cumulative distribution function (CDF): gives percentage of particles p_i whose radius r is less than a certain radius r_i. It requires p_0 = 0 and p_end = 1. The CDF's probabilities can be number, length, area or volume based. The cumulative number distribution function (CUMULATIVE_NUMBER_DISTRIBUTION) is used as PSD in MercuryDPM, as it can be interpreted as the probability that a particle's radius is less than r_i: CDF(r<r_i) = p_i.

Probability density function (PDF): p_i is the percentage of particles whose radius is between r_i-1 and r_i. It requires p_0=0 and sum(p_i)=1. The PDF's probabilities can also be number, length, area or volume based. PDF's are utilized to convert any type of PDF to a probability number density function (PROBABILITYDENSITY_NUMBER_DISTRIBUTION) which in a next step are converted to the default CUMULATIVE_NUMBER_DISTRIBUTION.

Sieve data: p_i is the percentage of particles whose radius is between r_i and r_i+1. It requires sum(p_i)=1. relation to PDF: p_i = p(PDF)_i-1. Sieve data is not yet used in this class.

Default distribution is the CUMULATIVE_NUMBER_DISTRIBUTION.

Member Enumeration Documentation

◆ TYPE

enum PSD::TYPE
strong

Enum class which stores the possible types of CDFs and PDFs. Particle size distributions can be represented by different probabilities based on number of particles, length of particles, surface area of particles or volume of particles.

Enumerator
CUMULATIVE_NUMBER_DISTRIBUTION 
CUMULATIVE_LENGTH_DISTRIBUTION 
CUMULATIVE_VOLUME_DISTRIBUTION 
CUMULATIVE_AREA_DISTRIBUTION 
PROBABILITYDENSITY_NUMBER_DISTRIBUTION 
PROBABILITYDENSITY_LENGTH_DISTRIBUTION 
PROBABILITYDENSITY_AREA_DISTRIBUTION 
PROBABILITYDENSITY_VOLUME_DISTRIBUTION 
72  {
73  CUMULATIVE_NUMBER_DISTRIBUTION,
74  CUMULATIVE_LENGTH_DISTRIBUTION,
75  CUMULATIVE_VOLUME_DISTRIBUTION,
76  CUMULATIVE_AREA_DISTRIBUTION,
77  PROBABILITYDENSITY_NUMBER_DISTRIBUTION,
78  PROBABILITYDENSITY_LENGTH_DISTRIBUTION,
79  PROBABILITYDENSITY_AREA_DISTRIBUTION,
80  PROBABILITYDENSITY_VOLUME_DISTRIBUTION
81  };

Constructor & Destructor Documentation

◆ PSD() [1/2]

PSD::PSD ( )

Constructor; sets everything to 0 or default.

Default constructor; sets every data member to 0 or default.

33 {
34  for (auto& momenta: momenta_)
35  momenta = 0;
36  index_ = 0;
37 }

References index_, and momenta_.

Referenced by copy().

◆ PSD() [2/2]

PSD::PSD ( const PSD other)

Copy constructor with deep copy.

Copy constructor

43 {
45  momenta_ = other.momenta_;
46  index_ = other.index_;
47 }

References index_, momenta_, and particleSizeDistribution_.

◆ ~PSD()

PSD::~PSD ( )
default

Destructor; default destructor.

Destructor. Since there are no pointers in this class, there is no need for any actions here.

Member Function Documentation

◆ computeCentralMomenta()

void PSD::computeCentralMomenta ( )

compute central momenta of the user defined PSD.

Compute the central momenta of the PSD from their respective raw momenta.

903 {
905  Mdouble mean = momenta_[1];
906  momenta_[5] += -5 * mean * momenta_[4] + 10 * mean * mean * momenta_[3]
907  - 10 * mean * mean * mean * momenta_[2] + 4 * mean * mean * mean * mean * mean;
908  momenta_[4] += -4 * mean * momenta_[3] + 6 * mean * mean * momenta_[2] - 3 * mean * mean * mean * mean;
909  momenta_[3] += -3 * mean * momenta_[2] + 2 * mean * mean * mean;
910  momenta_[2] += -mean * mean;
911 }

References computeRawMomenta(), and momenta_.

Referenced by computeStandardisedMomenta().

◆ computeRawMomenta()

void PSD::computeRawMomenta ( )

compute raw momenta of the user defined PSD.

Compute the raw momenta of the inserted PSD by converting it to a PDF, calculating the moments m_i according to \( m_i = \sum_{j=0}^n \) scaling it by the number of particles inserted into the simulation (moments are computed from a PROBABILITYDENSITY_NUMBER_DISTRIBUTION) and converting it back to a CDF. See Wikipedia for details.

883 {
885  for (size_t im = 0; im < momenta_.size(); ++im)
886  {
887  // prevent summing up of moments for each time step of the simulation
888  momenta_[im] = 0;
889  for (auto& it : particleSizeDistribution_)
890  {
891  momenta_[im] += std::pow(it.radius, im) * it.probability;
892  }
893  }
894  // zeroth-moment equals particle number for a PROBABILITYDENSITY_NUMBER_DISTRIBUTION
897 }

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityToCumulative(), getInsertedParticleNumber(), momenta_, and particleSizeDistribution_.

Referenced by computeCentralMomenta().

◆ computeStandardisedMomenta()

void PSD::computeStandardisedMomenta ( )

compute standardised momenta of the user defined PSD.

Compute the standardised momenta of the PSD from their respective central momenta.

917 {
919  Mdouble std = std::sqrt(momenta_[2]);
920  momenta_[3] /= std * std * std;
921  momenta_[4] /= std * std * std * std;
922  momenta_[5] /= std * std * std * std * std;
923 }

References computeCentralMomenta(), and momenta_.

◆ convertCumulativeToCumulativeNumberDistribution()

void PSD::convertCumulativeToCumulativeNumberDistribution ( TYPE  CDFType)

convert any other CDF to a CUMULATIVE_NUMBER_DISTRIBUTION.

converts any of the CDFTypes to a the default CUMULATIVE_NUMBER_DISTRIBUTION based on their TYPE.

Parameters
[in]PDFTypeType of the PDF: PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_AREA_DISTRIBUTION or PROBABILITYDENSITY_VOLUME_DISTRIBUTION, Where L = Length, A = Area and V = Volume.
Todo:
TP: NOT WORKING! If anyone knows how to do it feel free to add
617 {
618  Mdouble sum = 0;
619  switch (CDFType)
620  {
621  default:
622  logger(ERROR, "Wrong CDFType");
623  break;
625  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
626  {
627  // add conversion here
628 
629  // sum up probabilities
630  sum += it->probability;
631  }
632  // normalize
633  for (auto& p : particleSizeDistribution_)
634  {
635  p.probability /= sum;
636  }
637  break;
639  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
640  {
641  // add conversion here
642 
643  // sum up probabilities
644  sum += it->probability;
645  }
646  // normalize
647  for (auto& p : particleSizeDistribution_)
648  {
649  p.probability /= sum;
650  }
651  break;
653  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
654  {
655  // add conversion here
656 
657  // sum up probabilities
658  sum += it->probability;
659  }
660  // normalize
661  for (auto& p : particleSizeDistribution_)
662  {
663  p.probability /= sum;
664  }
665  break;
666  }
667 }

◆ convertCumulativeToProbabilityDensity()

void PSD::convertCumulativeToProbabilityDensity ( )

Converts a CDF to a PDF by derivation.

Convert any type of CDF to a PDF. Probabilities are derivated for each radius by substracting the CDF probabilities (i.e. pPDF_i = pCDF_i - pCDF_i-1).

510 {
511  // subtract cumulative probabilities
512  Mdouble probabilityOld = 0, probabilityCDF;
513  for (auto& it : particleSizeDistribution_)
514  {
515  probabilityCDF = it.probability;
516  it.probability -= probabilityOld;
517  probabilityOld = probabilityCDF;
518  }
519  // check whether probability density distribution is valid
521 }

References particleSizeDistribution_, and validateProbabilityDensityDistribution().

Referenced by computeRawMomenta(), getVolumeDx(), getVolumetricMeanRadius(), insertManuallyByVolume(), setPSDFromCSV(), and setPSDFromVector().

◆ convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution()

void PSD::convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution ( )

convert a PROBABILITYDENSITY_NUMBER_DISTRIBUTION to a PROBABILITYDENSITY_VOLUME_DISTRIBUTION.

converts a PROBABILITYDENSITY_NUMBER_DISTRIBUTION to a PROBABILITYDENSITY_VOLUME_DISTRIBUTION. Used for getVolumetricMeanRadius() and insertManuallyByVolume().

592 {
593  Mdouble sum = 0;
594  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
595  {
596  it->probability *=
597  0.25 * (square(it->radius) + square((it - 1)->radius)) * (it->radius + (it - 1)->radius);
598  // old conversion
599  //p.probability /= cubic(p.radius);
600  // sum up probabilities
601  sum += it->probability;
602  }
603  // normalize
604  for (auto& p : particleSizeDistribution_)
605  {
606  p.probability /= sum;
607  }
608 }

References particleSizeDistribution_, and mathsFunc::square().

Referenced by getVolumeDx(), getVolumetricMeanRadius(), and insertManuallyByVolume().

◆ convertProbabilityDensityToCumulative()

void PSD::convertProbabilityDensityToCumulative ( )

Converts a PDF to a CDF by integration.

Convert any type of PDF to a CDF. Probabilities are integrated for each radius by cumulatively summing them up (i.e. p_i = p_i + p_i-1).

495 {
496  // add up probabilities
497  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
498  {
499  it->probability = std::min(1.0, it->probability + (it - 1)->probability);
500  }
501  // check whether cumulative distribution is valid
503 }

References particleSizeDistribution_, and validateCumulativeDistribution().

Referenced by computeRawMomenta(), getVolumeDx(), insertManuallyByVolume(), setPSDFromCSV(), and setPSDFromVector().

◆ convertProbabilityDensityToProbabilityDensityNumberDistribution()

void PSD::convertProbabilityDensityToProbabilityDensityNumberDistribution ( TYPE  PDFType)

convert any PDF to a PROBABILITYDENSITY_NUMBER_DISTRIBUTION.

converts any of the PDFTypes to a PROBABILITYDENSITY_NUMBER_DISTRIBUTION based on their TYPE. This is a helper function which enables the convertProbabilityDensityToCumulative function to convert the psd into the default TYPE (CUMULATIVE_NUMBER_DISTRIBUTION).

Parameters
[in]PDFTypeType of the PDF: PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_AREA_DISTRIBUTION or PROBABILITYDENSITY_VOLUME_DISTRIBUTION, Where L = Length, A = Area and V = Volume.
530 {
531  Mdouble sum = 0;
532  switch (PDFType)
533  {
534  default:
535  logger(ERROR, "Wrong PDFType");
536  break;
538  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
539  {
540  it->probability /= 0.5 * (it->radius + (it - 1)->radius);
541  // old conversion
542  //p.probability /= p.radius;
543  // sum up probabilities
544  sum += it->probability;
545  }
546  // normalize
547  for (auto& p : particleSizeDistribution_)
548  {
549  p.probability /= sum;
550  }
551  break;
553  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
554  {
555  it->probability /=
556  1.0 / 3.0 * (square(it->radius) + it->radius * (it - 1)->radius + square((it - 1)->radius));
557  // old conversion
558  //p.probability /= square(p.radius);
559  // sum up probabilities
560  sum += it->probability;
561  }
562  // normalize
563  for (auto& p : particleSizeDistribution_)
564  {
565  p.probability /= sum;
566  }
567  break;
569  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
570  {
571  it->probability /=
572  0.25 * (square(it->radius) + square((it - 1)->radius)) * (it->radius + (it - 1)->radius);
573  // old conversion
574  //p.probability /= cubic(p.radius);
575  // sum up probabilities
576  sum += it->probability;
577  }
578  // normalize
579  for (auto& p : particleSizeDistribution_)
580  {
581  p.probability /= sum;
582  }
583  break;
584  }
585 }

References ERROR, logger, particleSizeDistribution_, PROBABILITYDENSITY_AREA_DISTRIBUTION, PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, and mathsFunc::square().

Referenced by insertManuallyByVolume(), setPSDFromCSV(), and setPSDFromVector().

◆ copy()

PSD * PSD::copy ( ) const

Creates a copy on the heap and returns a pointer.

Copy method; creates a copy on the heap and returns its pointer.

Returns
pointer to the copy on the heap
62 {
63 #ifdef DEBUG_CONSTRUCTOR
64  std::cout << "PSD::copy() const finished" << std::endl;
65 #endif
66  return new PSD(*this);
67 }

References PSD().

◆ cutHighSizeRatio()

void PSD::cutHighSizeRatio ( )

Check if the size ratio is too high and cut it.

Checks if the Size ratio of the PSD is too high and cuts the PSD at head and tail by 10 percent to avoid inaccurate results.

816 {
817  Mdouble SR = getSizeRatio();
818  if (SR > 100)
819  {
820  logger(WARN, "Size ratio > 100; Cutting the PSD to avoid inaccurate results");
821  cutoffCumulativeNumber(0.1, 0.9, 0.5);
822  }
823 }

References cutoffCumulativeNumber(), getSizeRatio(), logger, and WARN.

◆ cutoffAndSqueezeCumulative()

void PSD::cutoffAndSqueezeCumulative ( Mdouble  quantileMin,
Mdouble  quantileMax,
Mdouble  squeeze,
Mdouble  minPolydispersity = 0.1 
)

cutoff the PSD at given quantiles and make it less polydisperse by squeezing it.

Cuts off the CDF at given minimum and maximum quantiles, applies a minimum polydispersity at the base and squeezes the distribution to make it less polydisperse.

Parameters
[in]quantileMinundersize quantile to cut off the lower part of the CDF.
[in]quantileMaxoversize quantile to cut off the upper part of the CDF.
[in]squeezeapplies a squeezing factor ([0,1]) which determines the degree the PDF gets squeezed.
[in]minPolydispersityapplies a minimum of polydispersity ([0,1]) at the base of the CDF.
711 {
712  Mdouble r50 = 0.5 * PSD::getNumberDx(50);
713  // cut off
714  cutoffCumulativeNumber(quantileMin, quantileMax, minPolydispersity);
715  // squeeze psd
716  for (auto& p: particleSizeDistribution_)
717  {
718  p.radius = r50 + (p.radius - r50) * squeeze;
719  }
720 }

References cutoffCumulativeNumber(), getNumberDx(), and particleSizeDistribution_.

◆ cutoffCumulativeNumber()

void PSD::cutoffCumulativeNumber ( Mdouble  quantileMin,
Mdouble  quantileMax,
Mdouble  minPolydispersity = 0.1 
)

cutoff the PSD at given quantiles.

Cuts off the CDF at given minimum and maximum quantiles and applies a minimum polydispersity at the base.

Parameters
[in]quantileMinundersize quantile to cut off the lower part of the CDF.
[in]quantileMaxoversize quantile to cut off the upper part of the CDF.
[in]minPolydispersityApplies a minimum of polydispersity ([0,1]) at the base of the CDF.
676 {
677  Mdouble radiusMin = getRadiusByQuantile(quantileMin);
678  Mdouble radiusMax = getRadiusByQuantile(quantileMax);
679  // to get a minimum polydispersity at the base
680  Mdouble radiusMinCut = std::min(radiusMin * (1 + minPolydispersity), radiusMax);
681  // cut off min
682  while (particleSizeDistribution_.front().radius <= radiusMinCut)
683  {
685  }
686  particleSizeDistribution_.insert(particleSizeDistribution_.begin(), {radiusMinCut, quantileMin});
687  particleSizeDistribution_.insert(particleSizeDistribution_.begin(), {radiusMin, 0});
688  // cut off max
689  while (particleSizeDistribution_.back().radius >= radiusMax)
690  {
691  particleSizeDistribution_.pop_back();
692  }
693  Mdouble radiusMaxCut = std::max(radiusMax - (1 - minPolydispersity) * (radiusMax - particleSizeDistribution_.back()
694  .radius), radiusMin);
695  particleSizeDistribution_.push_back({radiusMaxCut, quantileMax});
696  particleSizeDistribution_.push_back({radiusMax, 1});
697  logger(INFO, "PSD was cut successfully and now has a size ratio of %", getSizeRatio());
698 }

References getRadiusByQuantile(), getSizeRatio(), INFO, logger, and particleSizeDistribution_.

Referenced by cutHighSizeRatio(), and cutoffAndSqueezeCumulative().

◆ decrementNParticlesPerClass()

void PSD::decrementNParticlesPerClass ( )

Decrement nParticlesPerClass_ counter.

Decrement the nParticlesPerClass_ counter of the last inserted class by 1. This is utilized when a particle is generated but failed to be placed into an insertionBoundary.

786 {
788 }

References index_, and nParticlesPerClass_.

◆ decrementVolumePerClass()

void PSD::decrementVolumePerClass ( Mdouble  volume)

Decrement volumePerClass_ counter.

Decrement the volumePerClass_ counter of the last inserted class by the volume of the particle. This is utilized when a particle is generated but failed to be placed into an insertionBoundary.

Parameters
[in]volumeVolume of the particle which failed to be inserted.
796 {
797  volumePerClass_[index_]-= volume;
798 }

References index_, and volumePerClass_.

◆ drawSample()

Mdouble PSD::drawSample ( )

Draw a sample radius from a CUMULATIVE_NUMBER_DISTRIBUTION.

Draws a sample probability of a real uniform distribution. A random number is generated in range [0,1] and by linear interpolation a suitable radius is returned. This function is only valid for CumulativeNumberDistributions as particles are drawn and not volumes, areas or lengths.

Returns
A Radius for a randomly assigned probability.
Todo:
TP: We should add a variable seed. Maybe make it definable by the user?
112 {
113  // draw a number between 0 and 1, uniformly distributed
115  static std::mt19937 gen(0);
116  static std::uniform_real_distribution<Mdouble> dist(0, 1);
117  Mdouble prob = dist(gen);
118  // find the interval [low,high] of the psd in which the number sits
119  auto high = std::lower_bound(particleSizeDistribution_.begin(), particleSizeDistribution_.end(), prob);
120  auto low = std::max(particleSizeDistribution_.begin(), high - 1);
121  Mdouble rMin = low->radius;
122  Mdouble rMax = high->radius;
123  Mdouble pMin = low->probability;
124  Mdouble pMax = high->probability;
125  // interpolate linearly between [low.radius,high.radius] (assumption: CDF is piecewise linear)
126  Mdouble a = (prob - pMin) / (pMax - pMin);
127  return a * rMin + (1 - a) * rMax;
128 }

References particleSizeDistribution_.

◆ getInsertedParticleNumber()

int PSD::getInsertedParticleNumber ( ) const

Get the number of particles already inserted into the simulation.

Gets the number of particles already inserted into the simulation by summing up the particles inserted in each class.

Returns
An integer containing the number of particles already inserted into the simulation.
867 {
868  int sum = 0;
869  for (auto& it: nParticlesPerClass_)
870  sum += it;
871  return sum;
872 }

References nParticlesPerClass_.

Referenced by computeRawMomenta().

◆ getMaxRadius()

Mdouble PSD::getMaxRadius ( ) const

Get largest radius of the PSD.

Gets the maximum radius of the PSD.

Returns
A double which corresponds to the maximum radius of the PSD.
839 {
840  return particleSizeDistribution_.back().radius;
841 }

References particleSizeDistribution_.

Referenced by getSizeRatio().

◆ getMinRadius()

Mdouble PSD::getMinRadius ( ) const

Get smallest radius of the PSD.

Gets the minimal radius of the PSD.

Returns
A double which corresponds to the minimal radius of the PSD.
830 {
831  return particleSizeDistribution_[0].radius;
832 }

References particleSizeDistribution_.

Referenced by getSizeRatio().

◆ getMomenta()

std::array< Mdouble, 6 > PSD::getMomenta ( ) const

get momenta of the user defined PSD.

Get momenta of the user defined particleSizeDistributionVector_ vector.

Returns
Array of momenta corresponding to the first six moments of the user defined PSD.
930 {
931  return momenta_;
932 }

References momenta_.

◆ getNumberDx()

Mdouble PSD::getNumberDx ( Mdouble  x) const

Calculate a certain diameter (e.g. D10, D50, D90, etc.) from a percentile x of the number based PSD.

Gets the diameter from a certain percentile of the number based PSD.

Returns
A double which is the diameter corresponding to a certain percentile.
Parameters
[in]xdouble which determines the obtained diameter as a percentile of the PSD.
728 {
729  return 2.0 * getRadiusByQuantile(x / 100);
730 }

References getRadiusByQuantile().

Referenced by cutoffAndSqueezeCumulative().

◆ getParticleSizeDistribution()

std::vector< PSD::RadiusAndProbability > PSD::getParticleSizeDistribution ( ) const

Get the PSD vector.

Gets the vector containing radii and probabilities of the PSD.

Returns
A vector containing radii and probabilities of the PSD.
848 {
850 }

References particleSizeDistribution_.

◆ getRadiusByQuantile()

Mdouble PSD::getRadiusByQuantile ( Mdouble  quantile) const

Calculate the quantile of the PSD.

gets the radius from a certain quantile of the PSD

Returns
A double which is the radius corresponding to a certain quantile of the PSD.
Parameters
[in]quantiledouble which determines the returned radius as a quantile of the PSD.
752 {
753  logger.assert_always(quantile <= 1 && quantile >= 0, "quantile is not between 0 and 1");
754  // find the quantile corresponding to the PSD
755  auto high = std::lower_bound(particleSizeDistribution_.begin(), particleSizeDistribution_.end(), quantile);
756  auto low = std::max(particleSizeDistribution_.begin(), high - 1);
757  if (high->probability == low->probability)
758  return high->radius;
759  else
760  return low->radius + (high->radius - low->radius) * (quantile - low->probability) /
761  (high->probability - low->probability);
762 }

References logger, and particleSizeDistribution_.

Referenced by cutoffCumulativeNumber(), getNumberDx(), and getVolumeDx().

◆ getSizeRatio()

Mdouble PSD::getSizeRatio ( ) const

get the size ratio (width) of the PSD.

Gets the size ratio (width) of the PSD defined by the ratio of minimum to maximum particle radius.

Returns
A double which corresponds to the size ratio of the PSD.
805 {
806  Mdouble rMin = getMaxRadius();
807  Mdouble rMax = getMinRadius();
808  return rMin / rMax;
809 }

References getMaxRadius(), and getMinRadius().

Referenced by cutHighSizeRatio(), cutoffCumulativeNumber(), and setPSDFromCSV().

◆ getVolumeDx()

Mdouble PSD::getVolumeDx ( Mdouble  x) const

Calculate a certain diameter (e.g. D10, D50, D90, etc.) from a percentile x of the volume based PSD.

Gets the diameter from a certain percentile of the volume based PSD.

Returns
A double which is the diameter corresponding to a certain percentile.
Parameters
[in]xdouble which determines the obtained diameter as a percentile of the PSD.

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution(), convertProbabilityDensityToCumulative(), and getRadiusByQuantile().

◆ getVolumetricMeanRadius()

Mdouble PSD::getVolumetricMeanRadius ( ) const

get a volumetric mean radius of the PSD.

Gets a radius such that a monodisperse system has the same number of particles as a polydisperse system. (i.e. mean += p_i * 0.5*(r_i^3 + r_i-1^3)

Returns
A double which corresponds to the volumetric mean radius.
770 {
771  PSD psd = *this;
774  Mdouble mean = 0;
775  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
776  mean += it->probability * 0.5 * (cubic(it->radius) + cubic((it - 1)->radius));
777  mean = pow(mean, 1. / 3.);
778  return mean;
779 }

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution(), mathsFunc::cubic(), and particleSizeDistribution_.

◆ insertManuallyByVolume()

Mdouble PSD::insertManuallyByVolume ( Mdouble  volume)

Draw sample radius manually per size class and check the volumeAllowed of each size class to insert the PSD as accurate as possible.

Draws a sample probability of a real uniform distribution within a given size class. A random number is generated in the size class range [r_i,r_i+1] and by linear interpolation a suitable radius is returned. from a CUMULATIVE_VOLUME_DISTRIBUTION the volumeAllowed of each class is checked to insert the PSD as accurate as possible. This function is only valid for cumulativeNumberDistributions as particles are drawn and not volumes, areas or lengths. Furthermore this insertion routine is most accurate for non-continuous particle insertion.

Parameters
[in]volumevolume of the geometry to be filled.
Returns
A Radius for a randomly assigned probability of a specific size class.
140 {
141  // initialize the particleNumberPerClass vector if empty.
142  if (nParticlesPerClass_.empty())
144  // initialize the particleNumberPerClass vector if empty.
145  if (volumePerClass_.empty())
147 
148  for (auto it = particleSizeDistribution_.end() - 1; it != particleSizeDistribution_.begin(); --it)
149  {
150  static std::mt19937 gen(0);
151  static std::uniform_real_distribution<Mdouble> dist(0, 1);
152  Mdouble prob = dist(gen);
153  // map the probability to the current class interval
154  prob = (it - 1)->probability + (it->probability - (it - 1)->probability) * prob;
155  // find the interval [low,high] of the psd in which the number sits
156  auto high = std::lower_bound(particleSizeDistribution_.begin(), particleSizeDistribution_.end(), prob);
157  auto low = std::max(particleSizeDistribution_.begin(), high - 1);
158  Mdouble rMin = low->radius;
159  Mdouble rMax = high->radius;
160  Mdouble pMin = low->probability;
161  Mdouble pMax = high->probability;
162  // interpolate linearly between [low.radius,high.radius] (assumption: CDF is piecewise linear)
163  index_ = std::distance(particleSizeDistribution_.begin(), high);
164  Mdouble a = (prob - pMin) / (pMax - pMin);
165  Mdouble rad = a * rMin + (1 - a) * rMax;
166  //\todo generalize this to work for non-spherical particles
167  // Compare inserted volume against volumeAllowed (only works for spherical at the moment)
168  Mdouble radVolume = 4.0 / 3.0 * constants::pi * rad * rad * rad;
171  Mdouble volumeAllowed = particleSizeDistribution_[index_].probability * volume;
174  if (volumePerClass_[index_] > volumeAllowed)
175  {
176  continue;
177  }
178  else if (radVolume + volumePerClass_[index_] > volumeAllowed)
179  {
180  Mdouble differenceVolumeLow = -(volumePerClass_[index_] - volumeAllowed);
181  Mdouble differenceVolumeHigh = radVolume + volumePerClass_[index_] - volumeAllowed;
182  Mdouble volumeRatio = differenceVolumeLow / differenceVolumeHigh;
183  // If volumeRatio > 1 it will insert anyways, because it should be no problem for the distribution
184  prob = dist(gen);
185  if (prob <= volumeRatio)
186  {
188  volumePerClass_[index_] += radVolume;
189  return rad;
190  }
191  else
192  {
193  continue;
194  }
195  }
196  else
197  {
199  volumePerClass_[index_] += radVolume;
200  return rad;
201  }
202 
203  }
204  return 0;
205 }

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution(), convertProbabilityDensityToCumulative(), convertProbabilityDensityToProbabilityDensityNumberDistribution(), index_, nParticlesPerClass_, particleSizeDistribution_, constants::pi, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, and volumePerClass_.

◆ printPSD()

void PSD::printPSD ( ) const

Prints radii and probabilities of the PSD vector.

Prints the radii [m] and probabilities [%] of the psd vector. It currently supports nanometers, micrometers, milimeters and meters as size units.

74 {
75  if (particleSizeDistribution_.front().radius > 1e-1)
76  {
77  for (const auto p : particleSizeDistribution_)
78  {
79  logger(INFO, "%m\t %\%", p.radius, p.probability * 100);
80  }
81  }
82  else if (particleSizeDistribution_.front().radius > 1e-4)
83  {
84  for (const auto p : particleSizeDistribution_)
85  {
86  logger(INFO, "%mm\t %\%", p.radius * 1000, p.probability * 100);
87  }
88  }
89  else if (particleSizeDistribution_.front().radius > 1e-7)
90  {
91  for (const auto p : particleSizeDistribution_)
92  {
93  logger(INFO, "%um\t %\%", p.radius * 1e6, p.probability * 100);
94  }
95  }
96  else
97  {
98  for (const auto p : particleSizeDistribution_)
99  {
100  logger(INFO, "%nm\t %\%", p.radius * 1e9, p.probability * 100);
101  }
102  }
103 }

References INFO, logger, and particleSizeDistribution_.

◆ setDistributionNormal()

void PSD::setDistributionNormal ( Mdouble  mean,
Mdouble  standardDeviation,
int  numberOfBins 
)

create a PSD vector for a normal distribution.

sets the particle size distribution to a discretised normal (gaussian) cumulative number distribution, which covers the range of 3 * standardDeviation (99,73% of all values covered).

Parameters
[in]meanDouble representing the mean of the particle size distribution
[in]standardDeviationDouble representing the standard deviation of the particle size distribution
[in]numberOfBinsInteger determining the number of bins (aka. particle size classes or resolution) of the particle size distribution
322 {
323  logger.assert_always(mean > 3 * standardDeviation,
324  "Reduce standardDeviation of your normal distribution to avoid negative radii; The mean "
325  "should be greater than 3*standardDeviation");
326  if (!particleSizeDistribution_.empty())
327  {
329  }
330  Mdouble radMin = mean - 3 * standardDeviation;
331  Mdouble radMax = mean + 3 * standardDeviation;
332  std::vector<Mdouble> radii = helpers::linspace(radMin, radMax, numberOfBins);
333  std::vector<Mdouble> probabilities;
334  for (int i = 0; i < radii.size(); i++)
335  {
336  Mdouble probability = 0.5 * (1 + erf((radii[i] - mean) / (sqrt(2) * standardDeviation)));
337  probabilities.push_back(probability);
338  }
339  for (int j = 0; j < radii.size(); j++)
340  {
341  particleSizeDistribution_.push_back({radii[j], probabilities[j]});
342  }
344 }

References constants::i, helpers::linspace(), logger, particleSizeDistribution_, and validateCumulativeDistribution().

Referenced by DistributionToPSDSelfTest::setupInitialConditions().

◆ setDistributionUniform()

void PSD::setDistributionUniform ( Mdouble  radMin,
Mdouble  radMax,
int  numberOfBins 
)

create a PSD vector for a uniform distribution.

sets the particle size distribution to a discretised uniform (linear) cumulative number distribution

Parameters
[in]radMinDouble representing The smallest particle radius of the particle size distribution
[in]radMaxDouble representing the biggest particle radius of the particle size distribution
[in]numberOfBinsInteger determining the number of bins (aka. particle size classes or resolution) of the particle size distribution
298 {
299  if (!particleSizeDistribution_.empty())
300  {
302  }
303  std::vector<Mdouble> probabilities = helpers::linspace(0.0, 1.0, numberOfBins);
304  std::vector<Mdouble> radii = helpers::linspace(radMin, radMax, numberOfBins);
305  // combine radii and probabilities
306  for (int i = 0; i < radii.size(); ++i)
307  {
308  particleSizeDistribution_.push_back({radii[i], probabilities[i]});
309  }
311 }

References constants::i, helpers::linspace(), particleSizeDistribution_, and validateCumulativeDistribution().

Referenced by BoundariesSelfTest::BoundariesSelfTest(), FluxBoundarySelfTest::FluxBoundarySelfTest(), CubeInsertionBoundary::set(), InsertionBoundarySelfTest::setupInitialConditions(), FullRestartTest::setupInitialConditions(), Chute::setupInitialConditions(), ChuteWithHopper::setupInitialConditions(), and T_protectiveWall::T_protectiveWall().

◆ setParticleSizeDistribution()

void PSD::setParticleSizeDistribution ( std::vector< RadiusAndProbability particleSizeDistribution)

set the PSD by a suitable vector.

sets the PSD from a vector of the PSD::RadiusAndProbability class; used to read in PSDs from a restart file.

Parameters
[in]particleSizeDistributionvector containing the radii and probabilities specifying the PSD.
857 {
858  particleSizeDistribution_ = particleSizeDistribution;
859 }

References particleSizeDistribution_.

◆ setPSDFromCSV()

void PSD::setPSDFromCSV ( const std::string &  fileName,
TYPE  PSDType,
bool  headings = false,
Mdouble  unitScalingFactorRadii = 1.0 
)

read in the PSD vector with probabilities and radii saved in a .csv file.

creates the PSD vector from probabilities and radii saved in a .csv file. Radii should be located at the first column and probabilities at the second column of the file. The Type of PSD will be converted to the default cumulative number distribution function (CUMULATIVE_NUMBER_DISTRIBUTION) for further processing.

Parameters
[in]fileNameName of the .csv file containing the radii and probabilities of the PSD.
[in]PSDTypeType of the PSD: CUMULATIVE_VOLUME_DISTRIBUTION, CUMULATIVE_NUMBER_DISTRIBUTION, CUMULATIVE_LENGTH_DISTRIBUTION, CUMULATIVE_AREA_DISTRIBUTION, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, PROBABILITYDENSITY_NUMBER_DISTRIBUTION, PROBABILITYDENSITY_LENGTH_DISTRIBUTION PROBABILITYDENSITY_AREA_DISTRIBUTION.
[in]headingsIf TRUE the file is assumed to have headings and the first row will be skipped. If FALSE the file has no headings and the file will be read in as is. Default is FALSE.
[in]unitScalingFactorRadiiScaling factor of radii to match SI-units.
426 {
427  // read in csv file using the csvReader class
428  csvReader csv;
429  csv.setHeader(headings);
430  csv.read(fileName);
431  std::vector<Mdouble> radii = csv.getFirstColumn(unitScalingFactorRadii);
432  std::vector<Mdouble> probabilities = csv.getSecondColumn(1.0);
433  logger.assert_always(radii.size() == probabilities.size(), "The radii and probabilities vector have to be the "
434  "same size");
435  // combine radii and probabilities
436  for (int i = 0; i < radii.size(); ++i)
437  {
438  particleSizeDistribution_.push_back({radii[i], probabilities[i]});
439  }
440  logger.assert_always(!particleSizeDistribution_.empty(), "PSD cannot be empty");
441  switch (PSDType)
442  {
443  // default case is a CUMULATIVE_NUMBER_DISTRIBUTION
444  default:
446  break;
453  break;
459  break;
466  break;
470  break;
475  break;
480  break;
485  break;
486  }
487  logger(INFO, "A PSD with size ratio % is now ready to be set", getSizeRatio());
488 }

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityToCumulative(), convertProbabilityDensityToProbabilityDensityNumberDistribution(), CUMULATIVE_AREA_DISTRIBUTION, CUMULATIVE_LENGTH_DISTRIBUTION, CUMULATIVE_VOLUME_DISTRIBUTION, csvReader::getFirstColumn(), csvReader::getSecondColumn(), getSizeRatio(), constants::i, INFO, logger, particleSizeDistribution_, PROBABILITYDENSITY_AREA_DISTRIBUTION, PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_NUMBER_DISTRIBUTION, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, csvReader::read(), csvReader::setHeader(), validateCumulativeDistribution(), and validateProbabilityDensityDistribution().

Referenced by MultiplePSDSelfTest::setupInitialConditions(), PSDManualInsertionSelfTest::setupInitialConditions(), and PSDSelfTest::setupInitialConditions().

◆ setPSDFromVector()

MERCURY_DEPRECATED void PSD::setPSDFromVector ( std::vector< RadiusAndProbability psdVector,
TYPE  PSDType 
)

Deprecated version of reading in PSDs from a vector.

creates the psd vector from radii and probabilities filled in by hand. The Type of PSD will be converted to the default cumulative number distribution function (CUMULATIVE_NUMBER_DISTRIBUTION) for further processing.

Parameters
[in]psdVectorVector containing radii and probabilities ([0,1]).
[in]PSDTypeType of the PSD: CUMULATIVE_VOLUME_DISTRIBUTION, CUMULATIVE_NUMBER_DISTRIBUTION, CUMULATIVE_LENGTH_DISTRIBUTION, CUMULATIVE_AREA_DISTRIBUTION, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, PROBABILITYDENSITY_NUMBER_DISTRIBUTION, PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_AREA_DISTRIBUTION.
Deprecated:
This is the old way of inserting PSDs. In the future use setPSDFromCSV().
360 {
361  particleSizeDistribution_ = psdVector;
362  switch (PSDType)
363  {
364  // default case is a cumulative number distribution (CND)
365  default:
367  break;
374  break;
381  break;
387  break;
392  break;
396  break;
401  break;
406  break;
407  }
408 }

References convertCumulativeToProbabilityDensity(), convertProbabilityDensityToCumulative(), convertProbabilityDensityToProbabilityDensityNumberDistribution(), CUMULATIVE_AREA_DISTRIBUTION, CUMULATIVE_LENGTH_DISTRIBUTION, CUMULATIVE_VOLUME_DISTRIBUTION, particleSizeDistribution_, PROBABILITYDENSITY_AREA_DISTRIBUTION, PROBABILITYDENSITY_LENGTH_DISTRIBUTION, PROBABILITYDENSITY_NUMBER_DISTRIBUTION, PROBABILITYDENSITY_VOLUME_DISTRIBUTION, validateCumulativeDistribution(), and validateProbabilityDensityDistribution().

◆ validateCumulativeDistribution()

void PSD::validateCumulativeDistribution ( )

Validates if a CDF starts with zero and adds up to unity.

Validates if a distribution is cumulative by first checking if the psd vector is empty and that the scaling of probabilities is in the range [0,1]. Also it checks if each consecutive value is higher than the latter value (i.e. assuming piecewise linear CDF: p_i > p_i-1). Further it replaces probabilities if the CDF does not start with zero or does not end with unity (i.e. p_0=0 and p_end=1)

214 {
215  // ensure interval of probabilities to be [0,1]
216  for (auto p = 0; p < particleSizeDistribution_.size(); ++p)
217  particleSizeDistribution_[p].probability =
218  particleSizeDistribution_[p].probability / particleSizeDistribution_.back().probability;
219  // check whether the distribution is cumulative
220  for (auto it = particleSizeDistribution_.begin() + 1; it != particleSizeDistribution_.end(); ++it)
221  {
222  logger.assert_always(it->probability >= (it - 1)->probability, "psd is not cumulative", true);
223  }
224  // cdf needs to start with a probability of zero
225  if (particleSizeDistribution_[0].probability != 0)
226  {
227  logger(INFO, "adding a zero at the beginning of the psd");
229  {std::max(1e-3 * particleSizeDistribution_[0].radius,
230  2 * particleSizeDistribution_[0].radius -
231  particleSizeDistribution_[1].radius), 0});
232  }
233  // cdf needs to end with a probability of one
234  if (particleSizeDistribution_.back().probability < 1)
235  {
236  logger(INFO, "adding a one at the end of the psd");
237  particleSizeDistribution_.push_back({particleSizeDistribution_.back().radius, 1});
238  }
239  // cut off equal subsequent values on the end of the cdf to remove zero size classes.
240  for (auto i = particleSizeDistribution_.end() - 1; i != particleSizeDistribution_.begin(); i--)
241  {
242  if (i->probability == (i - 1)->probability)
243  {
245  .end(), i));
246  }
247  else
248  {
249  break;
250  }
251  }
252 }

References constants::i, INFO, logger, and particleSizeDistribution_.

Referenced by convertProbabilityDensityToCumulative(), setDistributionNormal(), setDistributionUniform(), setPSDFromCSV(), and setPSDFromVector().

◆ validateProbabilityDensityDistribution()

void PSD::validateProbabilityDensityDistribution ( )

Validates if the integral of the PDF equals to unity.

Validates if a PDF is valid by first checking if the PDF starts with zero (i.e. p_0=0). Further it checks if the integral is equal to unity (i.e. assuming piecewise constant PDF: sum(p_i)=1).

259 {
260  Mdouble sum = 0;
261  // check if the psd starts with zero probability (i.e. p_0=0)
262  if (particleSizeDistribution_[0].probability != 0)
263  {
264  logger(INFO, "adding a zero at the beginning of PDF");
266  {std::max(1e-3 * particleSizeDistribution_[0].radius,
267  2 * particleSizeDistribution_[0].radius -
268  particleSizeDistribution_[1].radius), 0});
269  }
270  // sum up probabilities to check if it equals to unity (sum(p_i)=1)
271  for (auto& it : particleSizeDistribution_)
272  {
273  sum += it.probability;
274  }
275  // ensure interval of probabilities to be [0,1] by normalization.
276  for (auto& it : particleSizeDistribution_)
277  it.probability /= sum;
278  // if the sum of probabilities is not equal to unity, sum up the normalized probabilities again
279  if (sum - 1 > 1e-6)
280  {
281  sum = 0;
282  for (auto& it: particleSizeDistribution_)
283  {
284  sum += it.probability;
285  }
286  }
287  logger.assert_debug(sum - 1 < 1e-6, "PDF is not valid: Integral of PDF is not equal to unity");
288 }

References INFO, logger, and particleSizeDistribution_.

Referenced by convertCumulativeToProbabilityDensity(), setPSDFromCSV(), and setPSDFromVector().

Friends And Related Function Documentation

◆ operator< [1/2]

bool operator< ( const PSD::RadiusAndProbability l,
const PSD::RadiusAndProbability r 
)
friend

determines if a certain value of the PSD vector is lower than another one. Used for std::lower_bound()

Required to use std::lower_bound for finding when the probability is higher than a certain value.

Returns
TRUE if probability from a vector of type PSD::RadiusAndProbability is higher than a certain value from a vector of type PSD::RadiusAndProbability and FALSE in the opposite case.
940 {
941  return l.probability < r.probability;
942 }

◆ operator< [2/2]

bool operator< ( const PSD::RadiusAndProbability l,
Mdouble  r 
)
friend

determines if a certain value of the PSD vector is lower than a double.

required to use std::lower_bound for finding when the probability provided as a double is higher than a certain value.

Returns
TRUE if probability as double is higher than a certain value from a PSD::RadiusAndProbability vector and FALSE in the opposite case.
951 {
952  return l.probability < prob;
953 }

◆ operator<<

std::ostream& operator<< ( std::ostream &  os,
PSD::RadiusAndProbability p 
)
friend

Writes to output stream.

Writes to output stream. This function is used for restart files.

Returns
a reference to an output stream.
980 {
981  os << p.radius << ' ' << p.probability << ' ';
982  return os;
983 }

◆ operator==

Mdouble operator== ( PSD::RadiusAndProbability  l,
Mdouble  r 
)
friend

Determines if a certain value of the PSD vector is equal to a double.

Required to use std::distance to find the index of the PSD size class in which a particle has to be inserted

Returns
A double which determines the size class (radius) a particle will be inserted to.
960 {
961  return l.radius == r;
962 }

◆ operator>>

std::istream& operator>> ( std::istream &  is,
PSD::RadiusAndProbability p 
)
friend

Reads from input stream.

reads from input stream. This function is used for restart files.

Returns
a reference to an input stream.
969 {
970  is >> p.radius;
971  is >> p.probability;
972  return is;
973 }

Member Data Documentation

◆ index_

int PSD::index_
private

Integer which determines the class in which a particle has to be inserted for the manual insertion routine.

Referenced by decrementNParticlesPerClass(), decrementVolumePerClass(), insertManuallyByVolume(), and PSD().

◆ momenta_

std::array<Mdouble, 6> PSD::momenta_ {}
private

Array of doubles which stores the moments of a user defined discrete PROBABILITYDENSITY_NUMBER_DISTRIBUTION.

Todo:
TW can we make this a local variable instead of a class variable?

Referenced by computeCentralMomenta(), computeRawMomenta(), computeStandardisedMomenta(), getMomenta(), and PSD().

◆ nParticlesPerClass_

std::vector<int> PSD::nParticlesPerClass_
private

Vector of integers which represents the number of inserted particles in each size class. The classes in this vector are defined to contain the particles between size r_i and r_i-1. (e.g. size class 12 consists of particles between size class 12 and 11 of the PDF)

Referenced by decrementNParticlesPerClass(), getInsertedParticleNumber(), and insertManuallyByVolume().

◆ particleSizeDistribution_

◆ volumePerClass_

std::vector<Mdouble> PSD::volumePerClass_
private

Vector of doubles which stores the volume of inserted particles for each size class. This vector is used in the insertManuallyByVolume() function to check if the volumeAllowed per class is exceeded and thus no further particles should be added to a certain class. The classes in this vector are defined to contain the volume of particles between size r_i and r_i-1. (e.g. size class 12 consists of the particles' volume between size class 12 and 11 of the PDF)

Referenced by decrementVolumePerClass(), and insertManuallyByVolume().


The documentation for this class was generated from the following files:
PSD::getSizeRatio
Mdouble getSizeRatio() const
get the size ratio (width) of the PSD.
Definition: PSD.cc:804
mathsFunc::square
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
csvReader::read
void read(const std::string &filename)
Reads .csv files into Mercury.
Definition: csvReader.cc:40
PSD::TYPE::PROBABILITYDENSITY_AREA_DISTRIBUTION
@ PROBABILITYDENSITY_AREA_DISTRIBUTION
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:45
PSD::RadiusAndProbability::radius
Mdouble radius
Definition: PSD.h:89
PSD::volumePerClass_
std::vector< Mdouble > volumePerClass_
Definition: PSD.h:333
PSD::PSD
PSD()
Constructor; sets everything to 0 or default.
Definition: PSD.cc:32
PSD::momenta_
std::array< Mdouble, 6 > momenta_
Definition: PSD.h:317
PSD::convertProbabilityDensityToCumulative
void convertProbabilityDensityToCumulative()
Converts a PDF to a CDF by integration.
Definition: PSD.cc:494
PSD::validateProbabilityDensityDistribution
void validateProbabilityDensityDistribution()
Validates if the integral of the PDF equals to unity.
Definition: PSD.cc:258
logger
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
PSD::convertProbabilityDensityToProbabilityDensityNumberDistribution
void convertProbabilityDensityToProbabilityDensityNumberDistribution(TYPE PDFType)
convert any PDF to a PROBABILITYDENSITY_NUMBER_DISTRIBUTION.
Definition: PSD.cc:529
csvReader::setHeader
void setHeader(bool headings)
Set the boolean hasHeader_.
Definition: csvReader.cc:109
INFO
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
PSD::particleSizeDistribution_
std::vector< RadiusAndProbability > particleSizeDistribution_
Definition: PSD.h:311
PSD::validateCumulativeDistribution
void validateCumulativeDistribution()
Validates if a CDF starts with zero and adds up to unity.
Definition: PSD.cc:213
PSD::TYPE::PROBABILITYDENSITY_LENGTH_DISTRIBUTION
@ PROBABILITYDENSITY_LENGTH_DISTRIBUTION
PSD::RadiusAndProbability::probability
Mdouble probability
Definition: PSD.h:90
Mdouble
double Mdouble
Definition: GeneralDefine.h:34
PSD::convertCumulativeToProbabilityDensity
void convertCumulativeToProbabilityDensity()
Converts a CDF to a PDF by derivation.
Definition: PSD.cc:509
PSD::getInsertedParticleNumber
int getInsertedParticleNumber() const
Get the number of particles already inserted into the simulation.
Definition: PSD.cc:866
PSD::TYPE::CUMULATIVE_VOLUME_DISTRIBUTION
@ CUMULATIVE_VOLUME_DISTRIBUTION
ERROR
LL< Log::ERROR > ERROR
Error log level.
Definition: Logger.cc:53
PSD::index_
int index_
Definition: PSD.h:338
PSD::nParticlesPerClass_
std::vector< int > nParticlesPerClass_
Definition: PSD.h:324
PSD::convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution
void convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution()
convert a PROBABILITYDENSITY_NUMBER_DISTRIBUTION to a PROBABILITYDENSITY_VOLUME_DISTRIBUTION.
Definition: PSD.cc:591
helpers::linspace
std::vector< Mdouble > linspace(Mdouble a, Mdouble b, int N)
Definition: Helpers.cc:887
WARN
LL< Log::WARN > WARN
Warning log level.
Definition: Logger.cc:54
PSD::TYPE::PROBABILITYDENSITY_VOLUME_DISTRIBUTION
@ PROBABILITYDENSITY_VOLUME_DISTRIBUTION
csvReader::getSecondColumn
std::vector< Mdouble > getSecondColumn(Mdouble scalingFactor)
Get second column of a .csv file and return it as a double.
Definition: csvReader.cc:145
PSD::computeCentralMomenta
void computeCentralMomenta()
compute central momenta of the user defined PSD.
Definition: PSD.cc:902
PSD::getNumberDx
Mdouble getNumberDx(Mdouble x) const
Calculate a certain diameter (e.g. D10, D50, D90, etc.) from a percentile x of the number based PSD.
Definition: PSD.cc:727
PSD::TYPE::CUMULATIVE_LENGTH_DISTRIBUTION
@ CUMULATIVE_LENGTH_DISTRIBUTION
constants::i
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51
PSD::getMaxRadius
Mdouble getMaxRadius() const
Get largest radius of the PSD.
Definition: PSD.cc:838
PSD::cutoffCumulativeNumber
void cutoffCumulativeNumber(Mdouble quantileMin, Mdouble quantileMax, Mdouble minPolydispersity=0.1)
cutoff the PSD at given quantiles.
Definition: PSD.cc:675
PSD::getRadiusByQuantile
Mdouble getRadiusByQuantile(Mdouble quantile) const
Calculate the quantile of the PSD.
Definition: PSD.cc:751
PSD
Contains a vector with radii and probabilities of a user defined particle size distribution (PSD)
Definition: PSD.h:63
csvReader
Enables reading of .csv files into MercuryDPM.
Definition: csvReader.h:46
PSD::getMinRadius
Mdouble getMinRadius() const
Get smallest radius of the PSD.
Definition: PSD.cc:829
PSD::computeRawMomenta
void computeRawMomenta()
compute raw momenta of the user defined PSD.
Definition: PSD.cc:882
mathsFunc::cubic
T cubic(const T val)
calculates the cube of a number
Definition: ExtendedMath.h:115
PSD::TYPE::CUMULATIVE_AREA_DISTRIBUTION
@ CUMULATIVE_AREA_DISTRIBUTION
csvReader::getFirstColumn
std::vector< Mdouble > getFirstColumn(Mdouble scalingFactor)
Get first column of a .csv file and return it as a double.
Definition: csvReader.cc:128
PSD::TYPE::PROBABILITYDENSITY_NUMBER_DISTRIBUTION
@ PROBABILITYDENSITY_NUMBER_DISTRIBUTION