![]() |
ASPECT
|
Public Types | |
enum | KernelFunction { KernelFunction::gaussian, KernelFunction::triangular, KernelFunction::uniform, KernelFunction::cutoff_function_w1_dealii } |
Public Member Functions | |
ParticlePDF (const unsigned int granularity, const double bandwidth, const KernelFunction kernel_function) | |
ParticlePDF (const double bandwidth, const KernelFunction kernel_function) | |
void | fill_from_particle_range (const typename Particles::ParticleHandler< dim >::particle_iterator_range particle_range, const unsigned int n_particles_in_cell) |
void | insert_kernel_sum_from_particle_range (const Point< dim > reference_point, std::array< unsigned int, dim > table_index, const unsigned int n_particles_in_cell, const typename Particles::ParticleHandler< dim >::particle_iterator_range particle_range) |
void | add_value_to_function_table (const unsigned int x_index, const unsigned int y_index, const unsigned int z_index, const double input_value) |
void | add_value_to_function_table (const std::array< unsigned int, dim > &index_point, const double input_value) |
void | add_value_to_function_table (const double input_value, const types::particle_index reference_particle_id) |
double | evaluate_function_at_index (const unsigned int x_index, const unsigned int y_index, const unsigned int z_index) const |
void | compute_statistical_values () |
double | get_max () const |
double | get_min () const |
double | get_standard_deviation () const |
Private Member Functions | |
double | apply_selected_kernel_function (const double distance) const |
double | kernelfunction_uniform (const double distance) const |
double | kernelfunction_triangular (const double distance) const |
double | kernelfunction_gaussian (const double distance) const |
Private Attributes | |
Table< dim, double > | function_output_table |
double | bandwidth |
KernelFunction | kernel_function |
unsigned int | granularity |
double | max |
double | min |
double | standard_deviation |
double | mean |
small_vector< double > | function_output_vector |
bool | is_defined_per_particle |
types::particle_index | max_particle_index |
types::particle_index | min_particle_index |
This class computes the point-density function that describes the clustering of particles. ParticlePDF holds the values of the point-density function at every point at which the function is defined. The function is defined either at granularity^dim regular intervals throughout each cell, or at the position of each particle in the cell. ParticlePDF computes statistical values from the point-density function including the mean and standard deviation.
Definition at line 47 of file distribution.h.
|
strong |
The KernelFunctions
enum class is a data structure which contains the kernel functions available for use in the Kernel Density Estimator.
Enumerator | |
---|---|
gaussian | |
triangular | |
uniform | |
cutoff_function_w1_dealii |
Definition at line 56 of file distribution.h.
aspect::Particle::ParticlePDF< dim >::ParticlePDF | ( | const unsigned int | granularity, |
const double | bandwidth, | ||
const KernelFunction | kernel_function | ||
) |
This is the constructor for ParticlePDF. This constructor is called when creating a point-density function which is defined at regular intervals throughout a cell.
granularity | determines the number of reference points at which the point density function is computed. The point density function is only defined at these reference points. |
bandwidth | determines the bandwidth to be used in the kernel function. |
kernel_function | determines which kernel function to use when generating the point-density function. |
aspect::Particle::ParticlePDF< dim >::ParticlePDF | ( | const double | bandwidth, |
const KernelFunction | kernel_function | ||
) |
This is the constructor for ParticlePDF. This constructor is called when creating a point-density function which is defined at the positions of the particles within the cell, as opposed to at regular intervals based on a granularity value.
bandwidth | determines the bandwidth to be used in the kernel function. |
kernel_function | determines which kernel function to use when generating the point-density function. |
void aspect::Particle::ParticlePDF< dim >::fill_from_particle_range | ( | const typename Particles::ParticleHandler< dim >::particle_iterator_range | particle_range, |
const unsigned int | n_particles_in_cell | ||
) |
Fills the point-density function with values from the particles in the given cell.
particle_range | The particle_iterator_range to operate on. |
n_particles_in_cell | The number of particles belonging to the particle manager in question within the cell. |
void aspect::Particle::ParticlePDF< dim >::insert_kernel_sum_from_particle_range | ( | const Point< dim > | reference_point, |
std::array< unsigned int, dim > | table_index, | ||
const unsigned int | n_particles_in_cell, | ||
const typename Particles::ParticleHandler< dim >::particle_iterator_range | particle_range | ||
) |
This function is only called from fill_from_cell
. It iterates through every particle in the cell and sums the value of the kernel function between the reference point and the position of the cell.
reference_point | The point from which to get the value of the kernel function. |
table_index | The index in the PDF to insert the data into. |
n_particles_in_cell | The number of particles in the cell. |
particle_range | The particle_iterator_range to sum particles from. |
void aspect::Particle::ParticlePDF< dim >::add_value_to_function_table | ( | const unsigned int | x_index, |
const unsigned int | y_index, | ||
const unsigned int | z_index, | ||
const double | input_value | ||
) |
Inserts a value into the point-density function.
x_index | The x index in the point-density function to insert at. |
y_index | The y index in the point-density function to insert at. |
z_index | The z index in the point-density function to insert at. |
input_value | The value to insert into the point-density function. |
void aspect::Particle::ParticlePDF< dim >::add_value_to_function_table | ( | const std::array< unsigned int, dim > & | index_point, |
const double | input_value | ||
) |
Inserts a value into the point-density function.
index_point | The index in the point-density function to insert at. |
input_value | The value to insert into the point-density function. |
void aspect::Particle::ParticlePDF< dim >::add_value_to_function_table | ( | const double | input_value, |
const types::particle_index | reference_particle_id | ||
) |
Inserts a value into the point-density function. This version of the function is intended to be used when is_defined_per_particle
is true
.
input_value | The value to insert into the point-density function. |
reference_particle_id | the id of the particle from which the input value is taken in reference to. |
double aspect::Particle::ParticlePDF< dim >::evaluate_function_at_index | ( | const unsigned int | x_index, |
const unsigned int | y_index, | ||
const unsigned int | z_index | ||
) | const |
Returns the value of the function at the given index.
x_index | The x index in the point-density function to evaluate at. |
y_index | The y index in the point-density function to evaluate at. |
z_index | The z index in the point-density function to evaluate at. |
void aspect::Particle::ParticlePDF< dim >::compute_statistical_values | ( | ) |
Calculates the relevant statistics from the contents of the PDF.
double aspect::Particle::ParticlePDF< dim >::get_max | ( | ) | const |
Returns the maximum of the point-density function.
double aspect::Particle::ParticlePDF< dim >::get_min | ( | ) | const |
Returns the minimum of the point-density function.
double aspect::Particle::ParticlePDF< dim >::get_standard_deviation | ( | ) | const |
Returns the standard deviation of the point-density function.
|
private |
Returns the value of the selected kernel function.
distance | the distance to pass to the selected kernel function. |
coordinates | the coordinates representing the offset between the reference and sampled particle. This is needed because Functions::CutOffFunctionW1<dim> only takes a Point<dim> as input, not a double. |
|
private |
The Uniform kernel function returns a value of 1.0 as long as the distance is less than the KDE's bandwidth.
distance | the output of the kernel function depends on the distance between the reference point and the center of the kernel function. |
|
private |
The Triangular kernel function returns a value of 1.0 minus the distance variable.
distance | the output of the kernel function depends on the distance between the reference point and the center of the kernel function. |
|
private |
The gaussian function returns the value of a gaussian distribution at the specified distance.
distance | the output of the kernel function depends on the distance between the reference point and the center of the kernel function. |
|
private |
function_output_table
holds the output of the point-density function at every point where the function is defined, as long as the function has been generated at regular intervals throughout the cell (as opposed to being defined at the position of each particle.
Definition at line 191 of file distribution.h.
|
private |
The bandwidth
variable scales the point-density function. Choosing an appropriate bandwidth is important because a bandwidth value which is too low or too high can result in oversmoothing or undersmoothing of the point-density function. Oversmoothing or undersmoothing results in a function which represents the underlying data less accurately.
Definition at line 201 of file distribution.h.
|
private |
kernel_function
is an internal variable to keep track of which kernel function is being used by the ParticlePDF.
Definition at line 207 of file distribution.h.
|
private |
granularity
determines the number of inputs at which the point-density function is defined. The function is defined for granularity^dim inputs.
Definition at line 214 of file distribution.h.
|
private |
max
holds the maximum value of the point-density function after it has been computed.
Definition at line 219 of file distribution.h.
|
private |
min
holds the minimum value of the point-density function after it has been computed.
Definition at line 224 of file distribution.h.
|
private |
standard_deviation
holds the standard_deviation of the point-density function after it has been computed.
Definition at line 229 of file distribution.h.
|
private |
mean
holds the mean value of the point-density function after it has been computed.
Definition at line 234 of file distribution.h.
|
private |
function_output_vector
holds the output generated by summing the kernel functions from the positions of each particle in the cell, as opposed to summing the function at regular intervals throughout the cell.
Definition at line 241 of file distribution.h.
|
private |
is_defined_per_particle
keeps track of whether or not this point-density function is defined at regular intervals throughout the cell or if the function is defined at the position of each particle in the cell.
Definition at line 248 of file distribution.h.
|
private |
The index of the particle in the cell with the maximum point-density function value.
Definition at line 254 of file distribution.h.
|
private |
The index of the particle in the cell with the minimum point-density function value.
Definition at line 260 of file distribution.h.