Help entries
ComplexMixtures.CMTypes
— TypeUnion of types to define an isaprox
for testing.
ComplexMixtures.Box
— Typestruct Box{FloatVector, IntVector}
Structure that contains some data required to compute the linked cells.
sides::Any
nc::Any
l::Any
lcell::Int64
ComplexMixtures.ChemFile
— Typestruct ChemFile{T<:(AbstractVector)} <: Trajectory
Structure to contain a trajectory as read by Chemfiles.jl
filename::String
format::AbstractString
stream::Vector{Chemfiles.Trajectory}
nframes::Int64
sides::Vector{T} where T<:(AbstractVector)
solute::Selection
solvent::Selection
x_solute::Vector{T} where T<:(AbstractVector)
x_solvent::Vector{T} where T<:(AbstractVector)
natoms::Int64
ComplexMixtures.ChemFile
— MethodChemFile(filename::String, solute::Selection, solvent::Selection;format="" , T::Type = SVector{3,Float64})
Function open will set up the IO stream of the trajectory, fill up the number of frames field and additional parameters if required
The IO stream must be returned in a position such that the "nextframe!" function will be able to read the first frame of the trajectory
ComplexMixtures.CutoffDistances
— Typestruct CutoffDistances
Structure to contain a list of all the distances smaller than the cutoff.
nd::Vector{Int64}
d::Vector{Float64}
iat::Vector{Int64}
jat::Vector{Int64}
imol::Vector{Int64}
jmol::Vector{Int64}
maxdim::Vector{Int64}
ComplexMixtures.Density
— Typemutable struct Density
Structure to contain the density values obtained from the calculation.
solute::Float64
Default: 0.0
solvent::Float64
Default: 0.0
solvent_bulk::Float64
Default: 0.0
ComplexMixtures.DminMol
— Typemutable struct DminMol
Auxiliary structure used to update counters.
d::Float64
jmol::Int64
iat::Int64
jat::Int64
ComplexMixtures.FrameData
— Typemutable struct FrameData{T<:Trajectory, V}
Structure to contain data needed to compute the mddf for a single frame
trajectory::Trajectory
volume_frame::ComplexMixtures.Volume
rdf_count_random_frame::Vector{Float64}
md_count_random_frame::Vector{Float64}
dc::ComplexMixtures.CutoffDistances
dmin_mol::Vector{ComplexMixtures.DminMol}
dref_mol::Vector{Float64}
x_solvent_random::Vector
lc_solvent::ComplexMixtures.LinkedCells
ComplexMixtures.LinkedCells
— TypeLinkedCells
firstatom
is a vector that contains for each index i3D = index3D(i,j,k)
of a box, which is the first atom of the list of atoms in that box.
nextatom
is a vector that contains the index of the next atom of that box given the index of the previous atom (starting with firstatom
).
ComplexMixtures.NamdDCD
— Typestruct NamdDCD{T<:(AbstractVector)} <: Trajectory
Structure to contain the data of a trajectory in NAMD/DCD format.
filename::String
stream::FortranFiles.FortranFile
nframes::Int64
sides::Vector{T} where T<:(AbstractVector)
solute::Selection
solvent::Selection
x_solute::Vector{T} where T<:(AbstractVector)
x_solvent::Vector{T} where T<:(AbstractVector)
sides_in_dcd::Bool
lastatom::Int64
sides_read::Vector{Float64}
x_read::Vector{Float32}
y_read::Vector{Float32}
z_read::Vector{Float32}
ComplexMixtures.NamdDCD
— MethodNamdDCD(filename::String, solute::Selection, solvent::Selection;T::Type = SVector{3,Float64})
This function initializes the structure above, returning the data and the vectors with appropriate lengths and, importantly, with the i/o stream OPENED, ready to read the first trajectory frame using the "nextframe" function.
ComplexMixtures.Options
— Typestruct Options
Structure that contains the detailed input options.
firstframe::Int64
Default: 1
lastframe::Int64
Default: -1
stride::Int64
Default: 1
periodic::Bool
Default: true
irefatom::Int64
Default: -1
n_random_samples::Int64
Default: 10
binstep::Float64
Default: 0.02
dbulk::Float64
Default: 10.0
cutoff::Float64
Default: 10.0
usecutoff::Bool
Default: false
lcell::Int64
Default: 2
sleep::Float64
Default: 0.01
GC::Bool
Default: true
GC_threshold::Float64
Default: 0.1
seed::Int64
Default: 321
StableRNG::Bool
Default: false
nthreads::Int64
Default: -1
silent::Bool
Default: false
ComplexMixtures.OutputFiles
— Typemutable struct OutputFiles
Structure to contain the names of the output files.
output::String
solute_atoms::String
solvent_atoms::String
ComplexMixtures.Overview
— Typemutable struct Overview
Structure that is used to dispatch the show of a overview.
R::Result
domain_molar_volume::Float64
Default: 0.0
density::ComplexMixtures.Density
Default: Density()
solvent_molar_volume::Float64
Default: 0.0
solvent_molar_volume_bulk::Float64
Default: 0.0
solute_molar_volume::Float64
Default: 0.0
ComplexMixtures.PDBTraj
— Typestruct PDBTraj{T<:(AbstractVector)} <: Trajectory
Structure to contain PDB trajectories. Frames must be separated by "END", and with periodic cell sizes in the "CRYST1" field.
This structure and functions can be used as a template to implement the reading of other trajectory formats.
filename::String
stream::IOStream
nframes::Int64
sides::Vector{T} where T<:(AbstractVector)
solute::Selection
solvent::Selection
x_solute::Vector{T} where T<:(AbstractVector)
x_solvent::Vector{T} where T<:(AbstractVector)
natoms::Int64
x_read::Matrix{Float64}
ComplexMixtures.PDBTraj
— MethodPDBTraj(pdbfile::String, solute::Selection, solvent::Selection;T::Type = SVector{3,Float64})
Function open will set up the IO stream of the trajectory, fill up the number of frames field and additional parameters if required
The IO stream must be returned in a position such that the "nextframe!" function will be able to read the first frame of the trajectory
ComplexMixtures.Result
— Typestruct Result
Structure to contain the results of the MDDF calculation.
nbins::Int64
dbulk::Float64
cutoff::Float64
d::Vector{Float64}
Default: zeros(nbins)
md_count::Vector{Float64}
Default: zeros(nbins)
md_count_random::Vector{Float64}
Default: zeros(nbins)
sum_md_count::Vector{Float64}
Default: zeros(nbins)
sum_md_count_random::Vector{Float64}
Default: zeros(nbins)
mddf::Vector{Float64}
Default: zeros(nbins)
kb::Vector{Float64}
Default: zeros(nbins)
autocorrelation::Bool
solvent::ComplexMixtures.SolSummary
solute::ComplexMixtures.SolSummary
solute_atom::Matrix{Float64}
Default: zeros(nbins, solute.natomspermol)
solvent_atom::Matrix{Float64}
Default: zeros(nbins, solvent.natomspermol)
rdf_count::Vector{Float64}
Default: zeros(nbins)
rdf_count_random::Vector{Float64}
Default: zeros(nbins)
sum_rdf_count::Vector{Float64}
Default: zeros(nbins)
sum_rdf_count_random::Vector{Float64}
Default: zeros(nbins)
rdf::Vector{Float64}
Default: zeros(nbins)
kb_rdf::Vector{Float64}
Default: zeros(nbins)
density::ComplexMixtures.Density
Default: Density()
volume::ComplexMixtures.Volume
Default: Volume(nbins)
options::Options
irefatom::Int64
lastframe_read::Int64
nframes_read::Int64
files::Vector{String}
weights::Vector{Float64}
ComplexMixtures.Samples
— Typestruct Samples
Simple structure to contain the number of samples of each type of calculation to compute final results
md::Float64
random::Int64
ComplexMixtures.Selection
— Typestruct Selection
Structure that contains the information about the solute and solvent molecules.
natoms::Int64
nmols::Int64
natomspermol::Int64
index::Vector{Int64}
imol::Vector{Int64}
names::Vector{String}
ComplexMixtures.SolSummary
— Typestruct SolSummary
Structures to contain the details of a solute or solvent to store in the results of the MDDF calculation.
natoms::Int64
nmols::Int64
natomspermol::Int64
ComplexMixtures.Trajectory
— TypeTrajectory
Trajectory data type. Default reading with the Chemfiles infrastructure, except for DCD and PDB trajectory files, if the "PDBTraj" option is provided. See memory issue (https://github.com/chemfiles/Chemfiles.jl/issues/44)
ComplexMixtures.Units
— Typestruct Units
Unit conversions.
mole::Any
Default: 6.022140857e23
Angs3tocm3::Any
Default: 1.0e24
Angs3toL::Any
Default: 1.0e27
Angs3tocm3permol::Any
Default: mole / Angs3tocm3
Angs3toLpermol::Any
Default: mole / Angs3toL
SitesperAngs3tomolperL::Any
Default: Angs3toL / mole
ComplexMixtures.Volume
— Typemutable struct Volume
Structures to contain the volumes obtained from calculations.
total::Float64
bulk::Float64
domain::Float64
shell::Vector{Float64}
Base.isapprox
— MethodBase.isapprox(r1::T, r2::T; debug=false) where T <: CMTypes
Function to test if two runs offered similar results. Mostly used in the package testing routines
Base.merge
— Methodmerge(r::Vector{ComplexMixutures.Result})
This function merges the results of MDDF calculations obtained by running the same analysis on multiple trajectories, or multiple parts of the same trajectory. It returns a Result structure of the same type, with all the functions and counters representing averages of the set provided weighted by the number of frames read in each Result set.
Base.write
— Methodwrite(R::ComplexMixtures.Result, filename::String, solute::Selection, solvent::Selection)
Function to write the final results to output files as simple tables that are human-readable and easy to analyze with other software
If the solute and solvent selections are provides, pass on the atom names.
Base.write
— Methodwrite(R::ComplexMixtures.Result, filename::String;
solute_names::Vector{String} = ["nothing"],
solvent_names::Vector{String} = ["nothing"])
Optional passing of atom names.
ComplexMixtures.VMDselect
— MethodVMDselect(inputfile::String, selection::String; vmd="vmd" )
Select atoms using vmd selection syntax, with vmd in background
Returns the list of index (one-based) and atom names
Function to return the selection from a input file (topology, coordinates, etc), by calling VMD in the background.
ComplexMixtures.center_to_origin!
— Methodcenter_to_origin!(x::Vector{T}, center::T) where T
Translates atoms of vectory in x array such that center is in the origin. (x
is a vector of vectors).
ComplexMixtures.centerofcoordinates
— Methodcenterofcoordinates(coor::AbstractVector{T}) where T
Computes the center of coordinates of a vector.
ComplexMixtures.contrib
— Methodcontrib(s::Selection, atom_contributions::Array{Float64}, selection)
Extract the contribution of a given atom type selection from the solute or solvent atomic contributions to the MDDF.
s
here is the solute or solvent selection (type ComplexMixtures.Selection
) atom_contributions
is the R.solute_atom
or R.solvent_atom
arrays of the Result
structure, and the last argument is the selection of atoms from the solute to be considered, given as a list of indexes, list of atom names, vector of PDBTools.Atom
s, or a PDBTools.Residue
.
ComplexMixtures.cutoffdcell!
— Methodcutoffdcell!(cutoff::Float64,
iat::Int, xat::T,
x_solvent::Vector{T},
lc_solvent::LinkedCells,
box::Box,
i::Int, j::Int, k::Int,
dc::CutoffDistances) where T
Function that computes all distance of a point "xat" to the atoms of the solvent found in the linked cell corresponding to indexes i, j, and k
Modifies the data of dc
ComplexMixtures.cutoffdcell_self!
— Methodcutoffdcell_self!(cutoff::Float64,
iat::Int, xat::T,
x_solvent::Vector{T},
lc_solvent::LinkedCells,
box::Box,
i::Int, j::Int, k::Int,
dc::CutoffDistances,
solvent::Selection,
imol::Int) where T
Function that computes all distance of a point xat
to the atoms of the solvent found in the linked cell corresponding to indexes i, j, and k.
Modifies the data of dc
ComplexMixtures.cutoffdistances!
— Methodcutoffdistances!(cutoff::Float64,
x_solute::AbstractVector{T},
x_solvent::Vector{T},
lc_solvent::LinkedCells,
box::Box,
dc::CutoffDistances) where T
This routine that returns a list of the distances between atoms that are smaller than a specified cutoff, for a given set of coordinates.
Returns nd
, the number of distances smaller than the cutoff, and modifies dc
.
ComplexMixtures.cutoffdistances_self!
— Methodcutoffdistances_self!(cutoff::Float64,
x_solute::AbstractVector{T},
x_solvent::Vector{T},
lc_solvent::LinkedCells,
box::Box,
dc::CutoffDistances,
solvent::Selection, imol::Int) where T
This routine that returns a list of the distances between atoms that are smaller than a specified cutoff, for a given set of coordinates.
Returns nd
, the number of distances smaller than the cutoff, and modifies dc
.
ComplexMixtures.distance
— Methoddistance(x,y,sides)
Functions to compute Euclidean distances between two 3-dimensional vectors, with periodic boundary conditions given by sides
.
ComplexMixtures.distance
— Methoddistance(x,y)
Functions to compute Euclidean distances between 2 3-dimensional vectors.
ComplexMixtures.eulermat
— Methodeulermat(beta, gamma, theta, deg::String)
This routine was added because it defines the rotation in the "human" way, an is thus used to set the position of the fixed molecules. deg
can only be "degree"
, in which case the angles with be considered in degrees. If no deg
argument is provided, radians are used.
That means: beta
is a counterclockwise rotation around x
axis. gamma
is a counterclockwise rotation around y
axis. theta
is a counterclockwise rotation around z
axis.
ComplexMixtures.finalresults!
— Methodfinalresults!(R::Result, options::Options, trajectory::Trajectory, samples::Samples)
Function that computes the final results of all the data computed by averaging according to the sampling of each type of data, and converts to common units.
Computes also the final distribution functions and KB integrals
This function modified the values contained in the R data structure
ComplexMixtures.gr
— Methodgr(R::Result) = gr(R.d,R.rdf_count,R.density.solvent_bulk,R.options.binstep)
If a Result structure is provided without further details, use the rdf count and the bulk solvent density.
ComplexMixtures.gr
— Methodgr(r::Vector{Float64}, count::Vector{Float64}, density::Float64, binstep::Float64)
Computes the radial distribution function from the count data and the density.
This is exactly a conventional g(r) if a single atom was chosen as the solute and solvent selections.
Returns both the g(r) and the kb(r)
ComplexMixtures.grid3D
— Methodgrid3D(solute,solute_atoms,mddf_result,output_file; dmin=1.5, ddax=5.0, step=0.5)
This function builds the grid of the 3D density function and fills an array of mutable structures of type Atom, containing the position of the atoms of grid, the closest atom to that position, and distance.
solute
is a ComplexMixtuers.Selection
, defining the solute. solute_atoms
is the corresponding vector of PDBTools.Atom
s, and mddf_result
is the result of a mddf_result
calculation with the correspondign solute.
dmin
and dmax
define the range of distance where the density grid will be built, and step
defines how fine the grid must be. Be aware that fine grids involve usually a very large (hundreds of thousands points).
All parameters can be provides as keyword parameters.
Example
julia> using ComplexMixtures, PDBTools
julia> pdb = readPDB("./system.pdb");
julia> R = ComplexMixtures.load("./results.json");
julia> protein = select(pdb,"protein");
julia> solute = ComplexMixtures.Selection(protein,nmols=1);
julia> grid = ComplexMixtures.grid3D(solute=solute, solute_atoms=protein, mddf_result=R, output_file="grid.pdb");
grid
will contain a vector of Atom
s with the information of the MDDF at each grid point, and the same data will be written in the grid.pdb
file. This PDB file can be opened in VMD, for example, and contain in the beta
field the contribution of each protein residue to the MDDF at each point in space relative to the protein, and in the occupancy
field the distance to the protein. Examples of how this information can be visualized are provided in the user guide of ComplexMixtures
.
ComplexMixtures.icell1D
— Methodicell1D(nc::AbstractVector, i, j, k)
Returns the index of the linked cell, in the 1D representation, from its i,j,k coordinates.
ComplexMixtures.icell3D
— Methodicell3D
Returns the i,j,k coordinates of the cell from the coordinates of an atom and box properties
ComplexMixtures.inbulk
— Methodinbulk(d, R::Result)
Function that returns if a distance is in the bulk region or not, according to the options
ComplexMixtures.increase_size!
— MethodIncreases size of arrays in dc structure if needed
ComplexMixtures.initcells!
— Methodinitcells!(x::AbstractVector{T}, box::Box, lc::LinkedCells) where T
Function that initializes the linked cells by computing to each cell each atom belongs and filling up the firstatom and nexatom arrays.
Modifies the data in the lc structure
ComplexMixtures.isautocorrelation
— MethodJust check if solute and solvent are the same
ComplexMixtures.itype
— MethodGiven the index of the atom in the vector of coordinates of the solute or the solvent, returns the type of the atom, that is, the index of this atom within the molecule (goes from 1 to natomspermol)
ComplexMixtures.load
— Methodload(filename::String)
Function to load the json saved results file into the Result
data structure.
ComplexMixtures.mddf
— Methodmddf(trajectory::Trajectory, options::Options)
Function that computes the minimum-distance distribution function, atomic contributions, and KB integrals, given the Trajectory
structure of the simulation and, optionally, parameters given as a second argument of the Options
type. This is the main function of the ComplexMixtures
package.
Examples
julia> trajectory = ComplexMixtures.Trajectory("./trajectory.dcd",solute,solvent);
julia> results = ComplexMixtures.mddf(trajectory);
or
julia> options = Options(lastframe=1000);
julia> results = ComplexMixtures.mddf(trajectory,options);
ComplexMixtures.mddf_frame!
— Methodmddf_frame!(iframe::Int, framedata::FrameData, options::Options, RNG, R::Result)
Computes the MDDF for a single frame, modifies the data in the R
(type Result
) structure.
ComplexMixtures.mddf_frame_self!
— Methodmddf_frame_self!(iframe::Int, framedata::FrameData, options::Options, RNG, R::Result)
Computes the MDDF for a single frame, modifies the data in the R
(Result
) structure.
ComplexMixtures.mddf_linkedcells
— Methodmddf_linkedcells(trajectory::Trajectory, options::Options,samples::Samples, RNG, mddf_compute!)
Computes the MDDF using linked cells, serial version.
ComplexMixtures.mddf_linkedcells_parallel
— Methodmddf_linkedcells_parallel(trajectory::Trajectory, options::Options,samples::Samples, RNG, mddf_compute!)
Computes the MDDF using linked cells, in parallel.
ComplexMixtures.minimumdistance
— MethodReturns the the minimum distance between the points considered.
ComplexMixtures.move!
— Methodmove!(x::AbstractVector{T}, newcm::AbstractVector,beta, gamma, theta) where T
Translates and rotates a molecule according to the desired input center of coordinates and Euler rotations modifyies the vector x.
ComplexMixtures.overview
— Methodoverview(R::Result)
Function that outputs the volumes and densities in the most natural units.
ComplexMixtures.partialsort_cutoff!
— Methodpartialsort_cutoff!(x,cutoff; by = x -> x)
Function that reorders x vector by putting in the first positions the elements with values smaller than cutoff
ComplexMixtures.random_move!
— Methodrandom_move!(x_ref::AbstractVector{T},
irefatom::Int,
sides::T,
x_new::AbstractVector{T}, RNG) where T
Function that generates a new random position for a molecule.
The new position is returned in x
, a previously allocated array.
x_solvent_random
might be a view of the array that contains all the solvent molecules.
ComplexMixtures.save
— Methodsave(R::Result, filename::String)
Function to write the result data structure to a json file.
ComplexMixtures.setbin
— Methodsetbin(d,step)
Function that sets to which histogram bin a data point pertains simple, but important to keep consistency over all calls.
ComplexMixtures.shellradius
— Methodshellradius(i,step)
Compute the point in which the radius comprises half of the volume of the shell.
ComplexMixtures.sphereradiusfromshellvolume
— Methodsphereradiusfromshellvolume(volume,step)
Computes the radius that corresponds to a spherical shell of a given volume.
ComplexMixtures.sphericalshellvolume
— Methodsphericalshellvolume(i,step)
Computes the volume of the spherical shell defined within [(i-1)step,istep].
ComplexMixtures.sum!
— Methodsum!(R1::Result, R2::Result)
Sum the counts of two Results structures, adding the result to the first structure as in R1 = R1 + R2.
ComplexMixtures.title
— MethodPrint some information about the run.
ComplexMixtures.update_counters_frame!
— MethodUpdate the data with the data accumulated in a frame
ComplexMixtures.updatecounters!
— MethodFunction that updates the counters in R
and returns n_dmin_in_bulk
given the output of cutoffdistances
.
If the solute and solvent selections are provided, update mdcount, rdfcount and the atom-specific counters
returns: n_dmin_in_bulk
: number of molecules with all the atoms in the bulk n_dref_in_bulk
: number of molecules with the reference atom in the bulk
ComplexMixtures.viewmol
— Methodviewmol(i::Int, x::Vector{T}, n::Int) where T
Returns a view of a coordinate vector corresponding to the atoms of a molecule with index i. n is the number of atoms of the molecule.
ComplexMixtures.which_types
— Methodwhich_types(s::Selection, indexes::Vector{Int})
Function that returns the list of the indexes of the types of the atoms in a selection. For example, if a selection corresponds to a solvent of water molecules: There are three types, 1, 2, and 3, corresponding to the three atoms of the water molecule. If the indexes provided are, for instance, 11, 12, and 13, corresponding to a water molecule, this function will return 1, 2 and 3.
This is used to get equivalent-atom contributions to the distribution functions. For example, the input indexes span all water molecules, the output of this function will be still the three indexes corresponding to the three types of atoms that exist in a water molecule.
It is not possible to compute the contribution of one individual water molecule if the distribution function was computed for all molecules. Thus, the necessity to identify the types of atoms involved in a selection.
ComplexMixtures.wrap!
— MethodFunctions that wrap the coordinates
They modify the coordinates of the input vector.
function wrap!(x::AbstractVector{T}, sides::T, center::T) where T <: AbstractVector
Wrap to a given center of coordinates
ComplexMixtures.wrap!
— Methodwrap!(x::AbstractVector{T}, sides::T) where T <: AbstractVector
Wrap to origin.
ComplexMixtures.wrap_cell
— Methodwrap_cell(nc::AbstractVector, i::Int, j::Int, k::Int)
Given the indexes of a cell, return the periodic cell which correspondst to it, if the cell is outside the main box.
ComplexMixtures.writexyz
— Methodwritexyz(x::Vector{T}, file::String) where T <: AbstractVector
Print test xyz file.
``