WillardChandler¶
Module: willard-chandler¶
- class pytim.willard_chandler.WillardChandler(universe, group=None, alpha=2.0, radii_dict=None, mesh=2.0, symmetry='spherical', cluster_cut=None, include_zero_radius=False, cluster_threshold_density=None, extra_cluster_groups=None, centered=False, warnings=False, autoassign=True, **kargs)[source]¶
Identifies the dividing surface using the Willard-Chandler method NOTE that this method does not identify surface atoms
(Willard, A. P.; Chandler, D. J. Phys. Chem. B 2010, 114, 1954–1958)
- Parameters:
universe (Object) – The MDAnalysis Universe, MDTraj trajectory or OpenMM Simulation objects.
group (AtomGroup) – An AtomGroup, or an array-like object with the indices of the atoms in the group. Will identify the interfacial molecules from this group
alpha (float) – The width of the Gaussian kernel
mesh (float) – The grid spacing for the density calculation
group – 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 GROMOS 43a1) will be used.
cluster_cut (float) – Cutoff used for neighbors or density-based cluster search (default: None disables the cluster analysis)
cluster_threshold_density (float) – Number density threshold for the density-based cluster search. ‘auto’ determines the threshold automatically. Default: None uses simple neighbors cluster search, if cluster_cut is not None
extra_cluster_groups (Object) – Additional groups, to allow for mixed interfaces
include_zero_radius (bool) – if false (default) exclude atoms with zero radius from the surface analysis (they are always included in the cluster search, if present in the relevant group) to avoid some artefacts.
centered (bool) – Center the
group
warnings (bool) – Print warnings
Example:
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import * >>> >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> >>> radii = pytim_data.vdwradii(G43A1_TOP) >>> >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, fast=True) >>> R, _, _, _ = pytim.utilities.fit_sphere(inter.triangulated_surface[0]) >>> print ("Radius={:.3f}".format(R)) Radius=19.970
- property layers¶
The method does not identify layers.
Example:
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import * >>> >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, fast=True) >>> inter.layers <AtomGroup with 0 atoms>
- property autoassign¶
(bool) assign layers every time a frame changes
- property include_zero_radius¶
(bool) include atoms with zero radius in the analysis (excluded by default)
- writecube(filename='pytim.cube', group=None, sequence=False)[source]¶
Write to cube files (sequences) the volumentric density and the atomic positions.
- Parameters:
filename (str) – the file name
sequence (bool) – if true writes a sequence of files adding the frame to the filename
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import MICELLE_PDB >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, mesh=2.0) >>> inter.writecube('dens.cube') # writes on dens.cube >>> inter.writecube('dens.cube',group=g) # writes also particles >>> inter.writecube('dens.cube',sequence=True) # dens.<frame>.cube
- writeobj(filename='pytim.obj', sequence=False)[source]¶
Write to wavefront obj files (sequences) the triangulated surface
- Parameters:
filename (str) – the file name
sequence (bool) – if true writes a sequence of files adding the frame to the filename
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import MICELLE_PDB >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, mesh=2.0) >>> inter.writeobj('surf.obj') # writes on surf.obj >>> inter.writeobj('surf.obj',sequence=True) # surf.<frame>.obj
- property analysis_group¶
(AtomGroup) the group, the surface of which should be computed
- is_buried(pos)¶
Checks wether an array of positions are located below the first interfacial layer
- 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.
- prepare_box()¶
Before the analysis, pack every molecule into the box. Keep the original positions for latter use.
- reset_labels()¶
Reset labels before interfacial analysis
- property surface_cluster_cut¶
(float) distance cut to calculate the surface clusters
- 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 thesymmetry
option is different fromplanar
, thecentered='origin'
option is equivalent tocentered='middle'
.
- class pytim.willard_chandler.Writevtk(interface)[source]¶
- density(filename='pytim_dens.vtk', sequence=False)[source]¶
Write to vtk files the volumetric density.
- Parameters:
filename (str) – the file name
sequence (bool) – if true writes a sequence of files adding the frame to the filename
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import MICELLE_PDB >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, mesh=2.0)
>>> inter.writevtk.density('dens.vtk') # writes on dens.vtk >>> inter.writevtk.density('dens.vtk',sequence=True) # dens.<n>.vtk
- particles(filename='pytim_part.vtk', group=None, sequence=False)[source]¶
Write to vtk files the particles in a group.
- Parameters:
filename (str) – the file name
sequence (bool) – if true writes a sequence of files adding the frame to the filename
group (AtomGroup) – if None, writes the whole universe
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import MICELLE_PDB >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, mesh=2.0)
>>> # writes on part.vtk >>> inter.writevtk.particles('part.vtk') >>> # writes on part.<frame>.vtk >>> inter.writevtk.particles('part.vtk',sequence=True)
- surface(filename='pytim_surf.vtk', sequence=False)[source]¶
Write to vtk files the triangulated surface.
- Parameters:
filename (str) – the file name
sequence (bool) – if true writes a sequence of files adding the frame to the filename
>>> import MDAnalysis as mda >>> import pytim >>> from pytim.datafiles import MICELLE_PDB >>> u = mda.Universe(MICELLE_PDB) >>> g = u.select_atoms('resname DPC') >>> inter= pytim.WillardChandler(u, group=g, alpha=3.0, mesh=2.0) >>> inter.writevtk.surface('surf.vtk') # writes on surf.vtk >>> inter.writevtk.surface('surf.vtk',sequence=True) # surf.<n>.vtk