ChaconTarazona¶
Module: chacon_tarazona¶
-
class
pytim.chacon_tarazona.
ChaconTarazona
(universe, alpha=2.0, tau=1.5, group=None, radii_dict=None, max_layers=1, normal='guess', molecular=True, info=True, mesh=None, centered=False, warnings=False, autoassign=True, **kargs)¶ Identifies the dividing surface using the Chacon-Tarazona method
(Chacón, E.; Tarazona, P. Phys. Rev. Lett. 91, 166103, 2003) (Tarazona, P.; Chacón, E. Phys. Rev. B 70, 235407, 2004)
Parameters: - universe (Universe) – The MDAnalysis universe
- alpha (float) – Molecular scale cutoff
- tau (float) – Particles within this distance form the surface will be added during the self-consistent procedure.
- molecular (bool) – Switches between search of interfacial molecules / atoms (default: True)
- group (AtomGroup) – Compute the density using this group
- radii_dict (dict) – Dictionary with the atomic radii of the elements in the group. If None is supplied, the default one (from GROMOS43a1) will be used.
Example:
>>> import MDAnalysis as mda >>> import numpy as np >>> import pytim >>> from pytim.datafiles import WATER_GRO >>> u = mda.Universe(WATER_GRO) >>> g = u.select_atoms('name OW') >>> interface = pytim.ChaconTarazona(u,alpha=2.,tau=1.5,group=g,info=False,molecular=False) >>> interface.writepdb('CT.pdb',centered=True) >>> print (repr(interface.layers)) array([[<AtomGroup with 175 atoms>], [<AtomGroup with 159 atoms>]], dtype=object)
-
alpha
¶ (float) real space cutoff
-
autoassign
¶ (bool) assign layers every time a frame changes
-
cluster_cut
¶ (real) cutoff for phase identification
-
cluster_threshold_density
¶ (float) threshold for the density-based filtering
-
extra_cluster_groups
¶ (ndarray) additional cluster groups
-
info
¶ (bool) print additional information
-
itim_group
¶ (AtomGroup) the group, the surface of which should be computed
-
label_planar_sides
()¶ Assign to all layers a label (the beta tempfactor) that can be used in pdb files. Additionally, set the new layers and sides.
-
max_layers
¶ (int) maximum number of layers to be identified
-
molecular
¶ (bool) whether to compute surface atoms or surface molecules
-
multiproc
¶ (bool) use parallel implementation
-
radii_dict
¶ (dict) custom atomic radii
-
surfaces
¶ Surfaces associated to the interface
-
writepdb
(filename='layers.pdb', centered='no', group='all', multiframe=True, tempfactors=None)¶ Write the frame to a pdb file, marking the atoms belonging to the layers with different beta factors.
Parameters: - filename (str) – the output file name
- centered (str) – ‘origin’, ‘middle’, or ‘no’
- group (AtomGroup) – if ‘all’ is passed, use universe
- multiframe (bool) – append to pdb file if True
- tempfactors (ndarray) – use this array as temp (beta) factors
- Example: save the positions (centering the interface in the cell)
- without appending
>>> import pytim >>> import pytim.datafiles >>> import MDAnalysis as mda >>> from pytim.datafiles import WATER_GRO >>> u = mda.Universe(WATER_GRO) >>> interface = pytim.ITIM(u) >>> interface.writepdb('layers.pdb',multiframe=False)
- Example: save the positions without centering the interface. This
- will not shift the atoms from the original position
(still, they will be put into the basic cell).
The
multiframe
option set toFalse
will overwrite the file
>>> interface.writepdb('layers.pdb',centered='no')
Note that if GITIM is used, and the
symmetry
option is different from'planar'
, thecentered='origin'
option is equivalent tocentered='middle'
.