FOX.ff.lj_param
A module for estimating Lennard-Jones parameters.
Examples
>>> import pandas as pd
>>> from FOX import MultiMolecule, example_xyz, estimate_lennard_jones
>>> xyz_file: str = example_xyz
>>> atom_subset = ['Cd', 'Se', 'O']
>>> mol = MultiMolecule.from_xyz(xyz_file)
>>> rdf: pd.DataFrame = mol.init_rdf(atom_subset=atom_subset)
>>> param: pd.DataFrame = estimate_lennard_jones(rdf)
>>> print(param)
sigma (Angstrom) epsilon (kj/mol)
Atom pairs
Cd Cd 3.95 2.097554
Cd Se 2.50 4.759017
Cd O 2.20 3.360966
Se Se 4.20 2.976106
Se O 3.65 0.992538
O O 2.15 6.676584
Index
|
Estimate the Lennard-Jones \(\sigma\) and \(\varepsilon\) parameters using an RDF. |
|
Convert a distribution function into a free energy function. |
API
- FOX.ff.lj_param.estimate_lj(rdf, temperature=298.15, sigma_estimate='base')[source]
Estimate the Lennard-Jones \(\sigma\) and \(\varepsilon\) parameters using an RDF.
Given a radius \(r\), the Lennard-Jones potential \(V_{LJ}(r)\) is defined as following:
\[V_{LJ}(r) = 4 \varepsilon \left( \left( \frac{\sigma}{r} \right )^{12} - \left( \frac{\sigma}{r} \right )^6 \right )\]The \(\sigma\) and \(\varepsilon\) parameters are estimated as following:
\(\sigma\): The radii at which the first inflection point or peak base occurs in rdf.
\(\varepsilon\): The minimum value in of the rdf ree energy multiplied by \(-1\).
All values are calculated per atom pair specified in rdf.
- Parameters
rdf (
pandas.DataFrame
) – A radial distribution function. The columns should consist of atom-pairs.temperature (
float
) – The temperature in Kelvin.sigma_estimate (
str
) – Whether \(\sigma\) should be estimated based on the base of the first peak or its inflection point. Accepted values are"base"
and"inflection"
, respectively.
- Returns
A Pandas DataFrame with two columns,
"sigma"
(Angstrom) and"epsilon"
(kcal/mol), holding the Lennard-Jones parameters. Atom-pairs from rdf are used as index.- Return type
See also
MultiMolecule.init_rdf()
Initialize the calculation of radial distribution functions (RDFs).
get_free_energy()
Convert a distribution function into a free energy function.
- FOX.ff.lj_param.get_free_energy(distribution, temperature=298.15, unit='kcal/mol', inf_replace=nan)[source]
Convert a distribution function into a free energy function.
Given a distribution function \(g(r)\), the free energy \(F(g(r))\) can be retrieved using a Boltzmann inversion:
\[F(g(r)) = -RT * \text{ln} (g(r))\]Two examples of valid distribution functions would be the radial- and angular distribution functions.
- Parameters
distribution (array-like) – A distribution function (e.g. an RDF) as an array-like object.
temperature (
float
) – The temperature in Kelvin.inf_replace (
float
, optional) – A value used for replacing all instances of infinity (np.inf
).unit (
str
) – The to-be returned unit. See scm.plams.Units for a comprehensive overview of all allowed values.
- Returns
An array-like object with a free-energy function (kj/mol) of distribution.
- Return type
See also
MultiMolecule.init_rdf()
Initialize the calculation of radial distribution functions (RDFs).
MultiMolecule.init_adf()
Initialize the calculation of distance-weighted angular distribution functions (ADFs).