4.6.1. Polymer analysis — MDAnalysis.analysis.polymer

Author

Richard J. Gowers

Year

2015, 2018

Copyright

GNU Public License v3

This module contains various commonly used tools in analysing polymers.

class MDAnalysis.analysis.polymer.PersistenceLength(atomgroups, **kwargs)[source]

Calculate the persistence length for polymer chains

The persistence length is the length at which two points on the polymer chain become decorrelated. This is determined by first measuring the autocorrelation (\(C(n)\)) of two bond vectors (\(\mathbf{a}_i, \mathbf{a}_{i + n}\)) separated by \(n\) bonds

\[C(n) = \langle \cos\theta_{i, i+n} \rangle = \langle \mathbf{a_i} \cdot \mathbf{a_{i+n}} \rangle\]

An exponential decay is then fitted to this, which yields the persistence length

\[C(n) \approx \exp\left( - \frac{n \bar{l_B}}{l_P} \right)\]

where \(\bar{l_B}\) is the average bond length, and \(l_P\) is the persistence length which is fitted

Parameters
  • atomgroups (iterable) – List of AtomGroups. Each should represent a single polymer chain, ordered in the correct order.

  • verbose (bool, optional) – Show detailed progress of the calculation if set to True.

results.bond_autocorrelation

the measured bond autocorrelation

Type

numpy.ndarray

results.lb

the average bond length

New in version 2.0.0.

Type

float

lb

Alias to the results.lb.

Deprecated since version 2.0.0: Will be removed in MDAnalysis 3.0.0. Please use results.lb instead.

Type

float

results.x

length of the decorrelation predicted by lp

New in version 2.0.0.

Type

numpy.ndarray

results.lp

calculated persistence length

New in version 2.0.0.

Type

float

lp

Alias to the results.lp.

Deprecated since version 2.0.0: Will be removed in MDAnalysis 3.0.0. Please use results.lp instead.

Type

float

results.fit

the modelled backbone decorrelation predicted by lp

New in version 2.0.0.

Type

numpy.ndarray

fit

Alias to the results.fit.

Deprecated since version 2.0.0: Will be removed in MDAnalysis 3.0.0. Please use results.fit instead.

Type

float

See also

sort_backbone()

for producing the sorted AtomGroup required for input.

Example

from MDAnalysis.tests.datafiles import TRZ_psf, TRZ
import MDAnalysis as mda
from MDAnalysis.analysis import polymer
u = mda.Universe(TRZ_psf, TRZ)

# this system is a pure polymer melt of polyamide,
# so we can select the chains by using the .fragments attribute
chains = u.atoms.fragments

# select only the backbone atoms for each chain
backbones = [chain.select_atoms('not name O* H*') for chain in chains]

# sort the chains, removing any non-backbone atoms
sorted_backbones = [polymer.sort_backbone(bb) for bb in backbones]
persistence_length = polymer.PersistenceLength(sorted_backbones)

# Run the analysis, this will average over all polymer chains
# and all timesteps in trajectory
persistence_length = persistence_length.run()
print(f'The persistence length is: {persistence_length.results.lp}')

# always check the visualisation of this:
persistence_length.plot()

New in version 0.13.0.

Changed in version 0.20.0: The run method now automatically performs the exponential fit

Changed in version 1.0.0: Deprecated PersistenceLength.perform_fit() has now been removed.

Changed in version 2.0.0: Former results are now stored as results.bond_autocorrelation. lb, lp, fit are now stored in a MDAnalysis.analysis.base.Results instance.

plot(ax=None)[source]

Visualize the results and fit

Parameters

ax (matplotlib.Axes, optional) – if provided, the graph is plotted on this axis

Returns

ax

Return type

the axis that the graph was plotted on

MDAnalysis.analysis.polymer.fit_exponential_decay(x, y)[source]

Fit a function to an exponential decay

\[y = \exp\left(- \frac{x}{a}\right)\]
Parameters
  • x (array_like) – The two arrays of data

  • y (array_like) – The two arrays of data

Returns

a – The coefficient a for this decay

Return type

float

Notes

This function assumes that data starts at 1.0 and decays to 0.0

MDAnalysis.analysis.polymer.sort_backbone(backbone)[source]

Rearrange a linear AtomGroup into backbone order

Requires that the backbone has bond information, and that only backbone atoms are provided (ie no side chains or hydrogens).

Parameters

backbone (AtomGroup) – the backbone atoms, not necessarily in order

Returns

sorted_backbone – backbone in order, so sorted_backbone[i] is bonded to sorted_backbone[i - 1] and sorted_backbone[i + 1]

Return type

AtomGroup

New in version 0.20.0.