Package potentials_from_particle_insertion

Description

Python package with a set of functions for calculating the g(r) using the Widom test-particle insertion method, and/or for fitting the g(r) to solve for the pair-potential using particle coordinates from e.g. confocal microscopy or cryo-TEM.

The codes for this package were developed as part of the work in reference [1], where the iterative use of test-particle insertion to solve for the pair potential was based on [2] with correction for periodic or finite boundary conditions based on equations described in [3,4]

References

[1] Bransen, M. (2024). Measuring interactions between colloidal (nano)particles. PhD thesis, Utrecht University.

[2] Stones, A. E., Dullens, R. P. A., & Aarts, D. G. A. L. (2019). Model-Free Measurement of the Pair Potential in Colloidal Fluids Using Optical Microscopy. Physical Review Letters, 123(9), 098002. https://doi.org/10.1103/PhysRevLett.123.098002

[3] Markus Seserno (2014). How to calculate a three-dimensional g(r) under periodic boundary conditions. https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf

[4] Kopera, B. A. F., & Retsch, M. (2018). Computing the 3D Radial Distribution Function from Particle Positions: An Advanced Analytic Approach. Analytical Chemistry, 90(23), 13909–13914. https://doi.org/10.1021/acs.analchem.8b03157

License

MIT license

Copyright (c) 2024 Maarten Bransen

Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

Sub-modules

potentials_from_particle_insertion.multicomponent

This file contains codes for multicomponent compatible versions of functions, i.e. for more than one 'type' of particle …

Functions

def load_TPI_results(filename)

loads text file containing TPI results as written by save_tpi_results()

Parameters

filename : str
name of the TPI result file to load

Returns

rmin : float
lower bound on interparticle distance
rmax : float
upper bound on interparticle distance
dr : float
bin width of interparticle distance
bincentres : list of float
centre position of each interparticle distance bin
dhrdf : list of float
g(r) values from DH calculation for each bin
rdfs : list of list of float
g(r) values from TPI for each bin for each iteration
potentials : list of list of float
u(r) values for each bin for each iteration
counts : list of list of float
bin counts for each bin for each iteration
errors : list of float
chi squared error for each iteration
def rdf_dist_hist_2d(coordinates, rmin=0, rmax=10, dr=None, handle_edge='rectangle', boundary=None, density=None, quiet=False, neighbors_upper_bound=None, workers=1)

calculates g(r) via a 'conventional' distance histogram method for a set of 2D coordinate sets. Provided for convenience.

Parameters

coordinates : numpy.array or list-like of numpy.array
list of sets of coordinates, where each item along the 0th dimension is a N*2 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one TEM image, a time step from a video, etc.), with each element of the array of form [y,x]. Each set of coordinates is not required to have the same number of particles but must have the same particle density. When no boundaries or a single boundary is given, it is assumed all coordinate sets share these same boundaries.
rmin : float, optional
lower bound for the pairwise distance, left edge of 0th bin. The default is 0.
rmax : float, optional
upper bound for the pairwise distance, right edge of last bin. The default is 10.
dr : float, optional
bin width for the pairwise distance bins. The default is (rmax-rmin)/20
handle_edge : one of ['none','rectangle','periodic rectangle','circle', (list of) callable], optional

specifies how to correct for edge effects in the radial distribution function. These edge effects occur due to particles closer than rmax to any of the boundaries missing part of their neighbour shells in the dataset. The following options for correcting for this are available:

  • 'none' or None or False (all equivalent): do not correct for edge effects at all. Note that in order to calculate the particle density, cuboidal boundaries are assumed even when boundary is not specified. This can be overridden by explicitely giving the particle density using the density parameter.
  • 'rectangle': correct for the missing volume in a square or rectangular boundary.
  • 'periodic rectangle': like 'rectangle', except in periodic boundary conditions (i.e. one side wraps around to the other), based on ref [1].
  • 'circle': correct for missing volume in a spherical boundary.
  • a custom callable function (or list thereof) can be given to correct for arbitrary and/or mixed boundary conditions. This function must take three arguments: a numpy array of bin edges, an N×2 numpy array of coordinates and a boundary as specified in boundary, and return an N × len(bin edges)-1 numpy array with a value between 0 and 1 specifying the fraction of the volume of each circular shell that is within the boundary.

The default is 'rectangle'.

boundary : array-like, optional

positions of the walls that define the bounding box of the coordinates. The format must be given as ((ymin,ymax),(xmin,xmax)) when handle_edge is 'none', 'rectangle' or 'periodic rectangle'. When handle_edge='circle' the coordinates of the origin and radius of the bouding circle must be given as (y,x,radius). For custom boundary handling, boundary can be any format as required by the edge edge correction function(s) passed to handle_edge.

If all coordinate sets share the same boundary, a single such boundary can be given, otherwise a list of such array-likes of the same length as coordinates can be given for specifying boundaries of each set in coordinates separately. The default None which determines the smallest set of boundaries encompassing a set of coordinates.

density : float or callable or list of callable, optional
number density of particles in the box to use for normalizing the values. The default is the average density based on coordinates and boundary. When using custom edge correction in handle_edge it is possible to instead pass a (list of) callable that takes the number of particles and boundary and returns the number density.
quiet : bool, optional
if True, no output is printed to the terminal by this function call. The default is False.
neighbors_upper_bound : int, optional
upper limit on the number of neighbors expected within rmax from a particle. Useful for reducing memory consumption in datasets with dimensions much larger than rmax. The default is None, which takes the total number of particles in the coordinate set as upper limit.
workers : int, optional
number of workers to use for parallel processing during the neighbour detection step. If -1 is given all CPU threads are used. Note: this is ignored when handle_edge='periodic rectangle'. The default is 1.

Returns

rvals : numpy.array
bin-edges of the radial distribution function.
bincounts : numpy.array
values for the bins of the radial distribution function

References

[1] Markus Seserno (2014). How to calculate a three-dimensional g(r) under periodic boundary conditions. https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf

def rdf_dist_hist_3d(coordinates, rmin=0, rmax=10, dr=None, handle_edge='cuboid', boundary=None, density=None, quiet=False, neighbors_upper_bound=None, workers=1)

calculates g(r) via a 'conventional' distance histogram method for a set of 3D coordinate sets. Provided for convenience.

Parameters

coordinates : numpy.array or list-like of numpy.array
list of sets of coordinates, where each item along the 0th dimension is a N*3 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one z-stack, a time step from a video, etc.), with each element of the array of form [z,y,x]. Each set of coordinates is not required to have the same number of particles but must have the same particle density. When no boundaries or a single boundary is given, it is assumed all coordinate sets share these same boundaries.
rmin : float, optional
lower bound for the pairwise distance, left edge of 0th bin. The default is 0.
rmax : float, optional
upper bound for the pairwise distance, right edge of last bin. The default is 10.
dr : float, optional
bin width for the pairwise distance bins. The default is (rmax-rmin)/20
handle_edge : one of ['none','cuboid','periodic cuboid','sphere', (list of) callable], optional

specifies how to correct for edge effects in the radial distribution function. These edge effects occur due to particles closer than rmax to any of the boundaries missing part of their neighbour shells in the dataset. The following options for correcting for this are available:

  • 'none' or None or False (all equivalent): do not correct for edge effects at all. Note that in order to calculate the particle density, cuboidal boundaries are assumed even when boundary is not specified. This can be overridden by explicitely giving the particle density using density.
  • 'cuboid': correct for the missing volume in cuboidal boundary conditions, e.g. a 3D rectangular box with right angles. Based on ref. [1]
  • 'periodic cuboid': like 'cuboid', except in periodic boundary conditions (i.e. one side wraps around to the other). Based on ref. [2].
  • 'sphere': correct for missing volume in spherical boundary conditions
  • a custom callable function (or list thereof) can be given to correct for arbitrary and/or mixed boundary conditions. This function must take three arguments: a numpy array of bin edges, an N×3 numpy array of coordinates and a boundary as specified in boundary, and return an N × len(bin edges)-1 numpy array with a value between 0 and 1 specifying the fraction of the volume of each spherical shell that is within the boundary.

The default is 'cuboid'.

boundary : array-like, optional

positions of the walls that define the bounding box of the coordinates. The format must be given as ((zmin,zmax),(ymin,ymax),(xmin,xmax)) when handle_edge is 'none', 'cuboid' or 'periodic cuboid'. When handle_edge='sphere' the coordinates of the origin and radius of the bouding sphere must be given as (z,y,x,radius).

If all coordinate sets share the same boundary a single such boundary can be given, otherwise a list of such array-likes of the same length as coordinates can be given for specifying boundaries of each set in coordinates separately. The default is None which determines the smallest set of boundaries encompassing a set of coordinates.

density : float or callable or list of callable, optional
number density of particles in the box to use for normalizing the values. The default is the average density based on coordinates and boundary. When using custom edge correction in handle_edge it is possible to instead pass a (list of) callable that takes the number of particles and boundary and returns the number density.
quiet : bool, optional
if True, no output is printed to the terminal by this function call. The default is False.
neighbors_upper_bound : int, optional
upper limit on the number of neighbors expected within rmax from a particle. Useful for reducing memory consumption in datasets with dimensions much larger than rmax. The default is None, which takes the total number of particles in the coordinate set as upper limit.
workers : int, optional
number of workers to use for parallel processing during the neighbour detection step. If -1 is given all CPU threads are used. Note: this is ignored when handle_edge='periodic cuboid'. The default is 1.

Returns

rvals : numpy.array
bin-edges of the radial distribution function.
bincounts : numpy.array
values for the bins of the radial distribution function

References

[1] Kopera, B. A. F., & Retsch, M. (2018). Computing the 3D Radial Distribution Function from Particle Positions: An Advanced Analytic Approach. Analytical Chemistry, 90(23), 13909–13914. https://doi.org/10.1021/acs.analchem.8b03157

[2] Markus Seserno (2014). How to calculate a three-dimensional g(r) under periodic boundary conditions. https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf

def rdf_insertion_binned_2d(coordinates, pairpotential, rmin=0, rmax=10, dr=None, handle_edge='rectangle', boundary=None, pairpotential_binedges=None, n_ins=1000, interpolate=True, avoid_boundary=False, avoid_coordinates=False, neighbors_upper_bound=None, workers=1, testparticle_func=None, **kwargs)

Calculate g(r) from insertion of test-particles into sets of existing 2D coordinates, averaged over bins of width dr, and based on the pairwise interaction potential u(r) (in units of kT).

Parameters

coordinates : numpy.array or list-like of numpy.array
list of sets of coordinates, where each item along the 0th dimension is a N*2 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one TEM image, a time step from a video, etc.), with each element of the array of form [y,x]. Each set of coordinates is not required to have the same number of particles but must have the same particle density. When no boundaries or a single boundary is given, it is assumed all coordinate sets share these same boundaries.
pairpotential : iterable
list of values for the pairwise interaction potential. Must have length of len(pairpotential_binedges)-1 and be in units of thermal energy kT
rmin : float, optional
lower bound for the pairwise distance, left edge of 0th bin. The default is 0.
rmax : float, optional
upper bound for the pairwise distance, right edge of last bin. The default is 10.
dr : float, optional
bin width for the pairwise distance bins. The default is (rmax-rmin)/20
handle_edge : one of ['none','rectangle','periodic rectangle','circle', (list of) callable], optional

specifies how to correct for edge effects in the radial distribution function. These edge effects occur due to particles closer than rmax to any of the boundaries missing part of their neighbour shells in the dataset. The following options for correcting for this are available:

  • 'rectangle': correct for the missing volume in a square or rectangular boundary.
  • 'periodic' rectangle': like 'rectangle', except in periodic boundary conditions (i.e. one side wraps around to the other), based on ref [1].
  • 'circle': correct for missing volume in a spherical boundary.
  • a custom callable function (or list thereof) can be given to correct for arbitrary and/or mixed boundary conditions. This function must take three arguments: a numpy array of bin edges, an N×2 numpy array of coordinates and a boundary as specified in boundary, and return an N × len(bin edges)-1 numpy array with a value between 0 and 1 specifying the fraction of the volume of each spherical shell that is within the boundary.

The default is 'rectangle'.

boundary : array-like or list of such, optional

positions of the walls that define the bounding box of the coordinates. The format must be given as ((ymin,ymax),(xmin,xmax)) when handle_edge is `'rectangle' or 'periodic rectangle'. When handle_edge='circle' the coordinates of the origin and radius of the bouding circle must be given as (y,x,radius). For custom boundary handling, boundary can be any format as required by the edge edge correction function(s) passed to handle_edge.

If all coordinate sets share the same boundary, a single such boundary can be given, otherwise a list of such array-likes of the same length as coordinates can be given for specifying boundaries of each set in coordinates separately. The default is None which determines the smallest set of boundaries encompassing a set of coordinates.

pairpotential_binedges : iterable, optional
bin edges corresponding to the values in pairpotential. The default is None, which uses the bins defined by rmin, rmax and dr.
n_ins : int, optional
the number of test-particles to insert into each item in coordinates. The default is 1000.
interpolate : bool, optional
whether to use linear interpolation for calculating the interaction of two particles using the values in pairpotential. The default is True. If False, the nearest bin value is used.
avoid_boundary : bool, optional
if True, all test-particles are inserted at least rmax away from any of the surfaces defined in boundary to avoid effects of the finite volume of the bounding box. The default is False, which uses an analytical correction factor for missing volume of test-particles near the boundaries. This parameter is ignored when using custom edge handling.
avoid_coordinates : bool, optional
whether to insert test-particles at least rmin away from the center of any of the 'real' coordinates in coordinates. The default is False.
neighbors_upper_bound : int, optional
upper limit on the number of neighbors expected within rmax from a particle. Useful for reducing memory consumption in datasets with dimensions much larger than rmax. The default is None, which takes the total number of particles in the coordinate set as upper limit.
workers : int, optional
number of workers to use for parallel processing during the neighbour detection step. If -1 is given all CPU threads are used. Note: this is ignored when handle_edge='periodic rectangle'. The default is 1.
testparticle_func : callable or list thereof
function that generates test-particle coordinates in the bounding box, required when using custom edge handling (and ignored otherwise). When avoid_coordinates=True this function must take 4 arguments: a boundary as given in boundary, a N×2 numpy.array of coordinates to avoid, a float specifing the radius around the coordinates to avoid, and an integer giving the number of coordinates to generate. Otherwise, it must only take the boundary and the number of coordinates to generate as arguments. It must return an n_ins × 2 numpy array of coordinates.

Returns

rvals : numpy.array
bin-edges of the radial distribution function.
pair_correlation : numpy.array
values for the rdf / pair correlation function in each bin
counter : numpy.array
number of pair counts that contributed to the (mean) values in each bin

References

[1] Markus Seserno (2014). How to calculate a three-dimensional g(r) under periodic boundary conditions. https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf

def rdf_insertion_binned_3d(coordinates, pairpotential, rmin=0, rmax=10, dr=None, handle_edge='cuboid', boundary=None, pairpotential_binedges=None, n_ins=1000, interpolate=True, avoid_boundary=False, avoid_coordinates=False, neighbors_upper_bound=None, workers=1, testparticle_func=None, **kwargs)

Calculate g(r) from insertion of test-particles into sets of existing 2D coordinates, averaged over bins of width dr, and based on the pairwise interaction potential u(r) (in units of kT).

Parameters

coordinates : numpy.array or list-like of numpy.array
list of sets of coordinates, where each item along the 0th dimension is a N*3 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one z-stack, a time step from a video, etc.), with each element of the array of form [z,y,x]. Each set of coordinates is not required to have the same number of particles but must have the same particle density. When no boundaries or a single boundary is given, it is assumed all coordinate sets share these same boundaries.
pairpotential : iterable
list of values for the pairwise interaction potential. Must have length of len(pairpotential_binedges)-1 and be in units of thermal energy kT
rmin : float, optional
lower bound for the pairwise distance, left edge of 0th bin. The default is 0.
rmax : float, optional
upper bound for the pairwise distance, right edge of last bin. The default is 10.
dr : float, optional
bin width for the pairwise distance bins. The default is (rmax-rmin)/20
handle_edge : one of ['none','cuboid','periodic cuboid','sphere', (list of) callable], optional

specifies how to correct for edge effects in the radial distribution function. These edge effects occur due to particles closer than rmax to any of the boundaries missing part of their neighbour shells in the dataset. The following options for correcting for this are available:

  • 'none' or None or False (all equivalent): do not correct for edge effects at all. Note that in order to calculate the particle density, cuboidal boundaries are assumed even when boundary is not specified. This can be overridden by explicitely giving the particle density using density.
  • 'cuboid': correct for the missing volume in cuboidal boundary conditions, e.g. a 3D rectangular box with right angles. Based on ref. [1]
  • 'periodic' cuboid': like 'cuboid', except in periodic boundary conditions (i.e. one side wraps around to the other). Based on ref. [2].
  • 'sphere': correct for missing volume in spherical boundary conditions
  • a custom callable function (or list thereof) can be given to correct for arbitrary and/or mixed boundary conditions. This function must take three arguments: a numpy array of bin edges, an N×3 numpy array of coordinates and a boundary as specified in boundary, and return an N × len(bin edges)-1 numpy array with a value between 0 and 1 specifying the fraction of the volume of each spherical shell that is within the boundary.

The default is 'cuboid'.

boundary : array-like, optional

positions of the walls that define the bounding box of the coordinates. The format must be given as ((zmin,zmax),(ymin,ymax),(xmin,xmax)) when handle_edge is 'none', 'cuboid' or 'periodic cuboid'. When handle_edge='sphere' the coordinates of the origin and radius of the bouding sphere must be given as (z,y,x,radius).

If all coordinate sets share the same boundary a single such boundary can be given, otherwise a list of such array-likes of the same length as coordinates can be given for specifying boundaries of each set in coordinates separately. The default is None which determines the smallest set of boundaries encompassing a set of coordinates.

pairpotential_binedges : iterable, optional
bin edges corresponding to the values in pairpotential. The default is None, which uses the bins defined by rmin, rmax and dr.
n_ins : int, optional
the number of test-particles to insert into each item in coordinates. The default is 1000.
interpolate : bool, optional
whether to use linear interpolation for calculating the interaction of two particles using the values in pairpotential. The default is True. If False, the nearest bin value is used.
avoid_boundary : bool, optional
if True, all test-particles are inserted at least rmax away from any of the surfaces defined in boundary to avoid effects of the finite volume of the bounding box. The default is False, which uses an analytical correction factor for missing volume of test-particles near the boundaries. This parameter is ignored when using custom edge handling.
avoid_coordinates : bool, optional
whether to insert test-particles at least rmin away from the center of any of the 'real' coordinates in coordinates. The default is False.
neighbors_upper_bound : int, optional
upper limit on the number of neighbors expected within rmax from a particle. Useful for reducing memory consumption in datasets with dimensions much larger than rmax. The default is None, which takes the total number of particles in the coordinate set as upper limit.
workers : int, optional
number of workers to use for parallel processing during the neighbour detection step. If -1 is given all CPU threads are used. Note: this is ignored when handle_edge='periodic rectangle'. The default is 1.
testparticle_func : callable or list thereof
function that generates test-particle coordinates in the bounding box, required when using custom edge handling (and ignored otherwise). When avoid_coordinates=True this function must take 4 arguments: a boundary as given in boundary, a N×3 numpy.array of coordinates to avoid, a float specifing the radius around the coordinates to avoid, and an integer giving the number of coordinates to generate. Otherwise, it must only take the boundary and the number of coordinates to generate as arguments. It must return an n_ins×3 numpy array of coordinates.

Returns

rvals : numpy.array
bin-edges of the radial distribution function.
pair_correlation : numpy.array
values for the rdf / pair correlation function in each bin
counter : numpy.array
number of pair counts that contributed to the (mean) values in each bin

References

[1] Kopera, B. A. F., & Retsch, M. (2018). Computing the 3D Radial Distribution Function from Particle Positions: An Advanced Analytic Approach. Analytical Chemistry, 90(23), 13909–13914. https://doi.org/10.1021/acs.analchem.8b03157

[2] Markus Seserno (2014). How to calculate a three-dimensional g(r) under periodic boundary conditions. https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf

def rdf_insertion_exact_3d(coordinates, pairpotential, rmax, dr, boundary, pairpotential_binedges=None, gen_prob_reps=1000, shell_prob_reps=10, interpolate=True, use_numba=True, rmin=0)

Warning

This function underwent very minimal testing and is computationally inefficient. Use at your own risk.

calculate g(r) from particle insertion method using particle coordinates and pairwise interaction potential u(r) (in units of kT). Inserts test- particles at a specific r for every real particle. Implementation based on ref [1]

Parameters

coordinates : list of numpy.array of floats of shape nt*n*3
The coordinates of particles in the system
pairpotential : list of floats
the value of the interparticle potential in each bin of 0 to rmax+dr in steps of dr. The value is assumed to correspond to the centrepoint of the bins.
rmax : float
cut-off value for the radius. Interactions beyond this range are considered negligable. Right edge of last bin in the pair potential
dr : float
bin width/step size of the interparticle distance used for g(r) and the pair potential.
boundary : tuple of floats of form ((zmin,zmax),(ymin,ymax),(xmin,xmax))
defines the boundaries of the box in which coordinates may exist.
gen_prob_reps : int, optional
number of trial particles used to evaluate the general probability of placing a particle at random coordinates in the box defined by boundary. The default is 1000.
shell_prob_reps : int, optional
number of trial particles per reference coordinate in coordinates used to evaluate the probability of placing a particle in each distance bin. The default is 10.
interpolate : bool, optional
If true, the pair potential is linearly interpolated between bin centres to calculate the interactions between all pairs of particles. If false, the value from the nearest bin is taken which is slower to compute. The default is True.

Returns

pair_correlation : list of float
the pair correlation / radial distribution functions evaluated in the bins whose edges are defined by numpy.arange(0,rmax+dr,dr)
counters : list of int
The number of trialparticles evaluated for each distance bin

References

[1] Stones, A. E., Dullens, R. P. A., & Aarts, D. G. A. L. (2018). Contact values of pair distribution functions in colloidal hard disks by test-particle insertion. The Journal of Chemical Physics, 148(24), 241102. https://doi.org/10.1063/1.5038668

def run_iteration(coordinates, pair_correlation_func, initial_guess=None, rmin=0, rmax=10, dr=None, convergence_tol=1e-05, max_iterations=100, zero_clip=1e-20, regulate=False, **kwargs)

Run the algorithm to solve for the pairwise potential that most accurately reproduces the radial distribution function using test-particle insertion, based on the ideas in ref. [1].

Parameters

coordinates : list-like of numpy.array
List of sets of coordinates, where each item along the 0th dimension is a n*3 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one z-stack, a time step from a video, etc.), with each element of the array of form [z,y,x] or [y,x] in case of 2D data. Each set of coordinates is not required to have the same number of particles but all stacks must share the same bounding box as given by boundary, and all coordinates must be within this bounding box.
pair_correlation_func : list of float
bin values for the true pair correlation function that the algorithm will try to match iteratively.
initial_guess : list of float, optional
Initial guess for the particle potential on the 0th iteration. The default is None which gives 0 in each bin.
rmin : float, optional
left edge of the smallest bin in interparticle distance r to consider. The default is 0.
rmax : float, optional
Right edge of the largest bin in interparticle distance r to consider. The default is 20.
dr : float, optional
Stepsize or bin width in interparticle distance r. The default is 0.5.
convergence_tol : float, optional
target value for χ², if it dips below this value the iteration is considered to be converged and ended. The default is 1e-5.
max_iterations : int, optional
Maximum number of iterations after which the algorithm is ended. The default is 100.
zero_clip : float, optional
values below the value of zero-clip are set to this value to avoid devision by zero errors. The default is 1e-20.
regulate : bool, optional
if True, use regularization to more gently nudge towards the input g(r) at the cost of slower convergence. Experimental option. The default is False.
**kwargs : key=value
Additional keyword arguments are passed on to rdf_insertion_binned_2d() or rdf_insertion_binned_3d()

Returns

χ² : list of float
summed squared error in the pair correlation function for each iteration
pairpotential : list of list of float
the values for the pair potential in each bin for each iteration
paircorrelation : list of list of float
the values for the pair correlation function from test-particle insertion for each iteration
counts : list of list of int
number of pair counts contributing to each bin in each iteration

References

[1] Stones, A. E., Dullens, R. P. A., & Aarts, D. G. A. L. (2019). Model- Free Measurement of the Pair Potential in Colloidal Fluids Using Optical Microscopy. Physical Review Letters, 123(9), 098002. https://doi.org/10.1103/PhysRevLett.123.098002

See Also

rdf_insertion_binned_2d()
2D routine for g(r) from test-particle insertion
rdf_insertion_binned_3d()
3D routine for g(r) from test-particle insertion
def run_iterator_fitfunction(coordinates, pair_correlation_func, potential_func, initial_guess=None, fit_bounds=None, rmin=0, rmax=20, dr=None, convergence_tol=1e-05, max_iterations=100, max_func_evals=100, zero_clip=1e-20, **kwargs)

Run the algorithm to solve for the pairwise potential that most accurately reproduces the radial distribution function using test-particle insertion, using a fit function of known form instead of a generalized binned potential based on scipy.optimize.curve_fit.

Parameters

coordinates : list-like of numpy.array
List of sets of coordinates, where each item along the 0th dimension is a n*3 numpy.array of particle coordinates, where each array is an independent set of coordinates (e.g. one z-stack, a time step from a video, etc.), with each element of the array of form [z,y,x] or [y,x] in case of 2D data. Each set of coordinates is not required to have the same number of particles but all stacks must share the same bounding box as given by boundary, and all coordinates must be within this bounding box.
pair_correlation_func : list of float
bin values for the true pair correlation function that the algorithm will try to match iteratively.
potential_func : callable
fit funtion to use, where the first argument is the interparticle distance r, and any subsequent arguments are optimized for.
initial_guess : list, optional
values for the function parameters to use for the first iteration. The default is None.
fit_bounds : 2-tuple of list_like, optional
Each element of the tuple must be either an array with the length equal to the number of parameters, or a scalar (in which case the bound is taken to be the same for all parameters). Use np.inf with an appropriate sign to disable bounds on all or some parameters.
rmin : float, optional
left edge of the smallest bin in interparticle distance r to consider. The default is 0.
rmax : float, optional
Right edge of the largest bin in interparticle distance r to consider. The default is 20.
dr : float, optional
Stepsize or bin width in interparticle distance r. The default is 0.5.
convergence_tol : float, optional
target value for χ², if it dips below this value the iteration is considered to be converged and ended. The default is 1e-5.
max_iterations : int, optional
Maximum number of iterations after which the algorithm is ended. The default is 100.
max_func_evals : int, optional
Maximum number of times the TPI function is evaluated. The default is 100.
zero_clip : float, optional
values below the value of zero-clip are set to this value to avoid devision by zero errors. The default is 1e-20.
**kwargs : key=value
Additional keyword arguments are passed on to rdf_insertion_binned_2d() or rdf_insertion_binned_3d()

Returns

χ² : list of float
summed squared error in the pair correlation function for each iteration
pairpotential : list of callable
potential functions giving the potential for each iteraction
fit_params : list of list of float
the values for the function parameters in each iteration
paircorrelation : list of list of float
the values for the pair correlation function from test-particle insertion for each iteration
counts : list of list of int
number of pair counts contributing to each bin in each iteration

See Also

rdf_insertion_binned_2d()
2D routine for g(r) from test-particle insertion
rdf_insertion_binned_3d()
3D routine for g(r) from test-particle insertion
def save_TPI_results(rmin, rmax, dr, bincentres, dhrdf, rdfs, potentials, counts, errors, filename='particle_insertion_result.txt')

utility function for saving results from TPI as txt file

Parameters

rmin : float
lower bound on interparticle distance
rmax : float
upper bound on interparticle distance
dr : float
bin width of interparticle distance
bincentres : list of float
centre position of each interparticle distance bin
dhrdf : list of float
g(r) values from DH calculation for each bin
rdfs : list of list of float
g(r) values from TPI for each bin for each iteration
potentials : list of list of float
u(r) values for each bin for each iteration
counts : list of list of float
bin counts for each bin for each iteration
errors : list of float
chi squared error for each iteration
filename : str, optional
filename to use. The default is 'particle_insertion_result.txt'.