Developer zone

Simulation

MolSimToolkit.SimulationType
Simulation(pdb_file::String, trajectory_file::String; first=1, last=nothing, step=1)
Simulation(atoms::AbstractVector{<:AtomType}, trajectory_file::String; first=1, last=nothing, step=1)

Creates a new Simulation object.

The first constructor creates a Simulation object from a PDB or mmCIF file and a trajectory file. It will use the PDBTools.Atom for the atom type, which will populate the atoms vector of the Simulation object. Currently, other atom types are supported, if the MolSimToolkit.atomic_mass(::AtomType) function is defined for the atom type.

With the second constructor, the atoms vector is passed as an argument. This is useful when the atoms are provided by a different source than the PDB file.

If first, last, and step are not specified, the Simulation will iterate over all frames in the file.

A Simulation object contains a trajectory file and a PDB data of the atoms. It can be iterated over to obtain the frames in the trajectory. The Simulation object is a mutable struct that contains the following data, that can be retrieved by the corresponding functions:

  • frame_range(::Simulation): the range of frames to be iterated over
  • frame_index(::Simulation): the index of the current frame in the trajectory
  • length(::Simulation): the number of frames to be iterated over in the trajectory file, considering the current range
  • raw_length(::Simulation): the number of frames in the trajectory file
  • atoms(::Simulation): the atoms in the simulation

The Simulation object can also be manipulated by the following functions:

  • close(::Simulation): closes the trajectory file
  • restart!(::Simulation): restarts the iteration over the trajectory file
  • first_frame!(::Simulation): restarts the iteration over the trajectory file and places the current frame at the first frame in the trajectory
  • current_frame(::Simulation): returns the current frame in the trajectory
  • next_frame!(::Simulation): reads the next frame in the trajectory file and returns it. Moves the current frame to the next one.
  • set_frame_range!(::Simulation; first, last, step): resets the range of frames to be iterated over.
  • get_frame(::Simulation, iframe): returns the frame at the given index in the trajectory.

One important feature of the Simulation object is that it can be iterated over, frame by frame.

The pairs iterator can also be used to iterate over the frames, returning a tuple with the frame index and the frame itself.

The enumerate iterator can also be used to iterate over the frames, returning a tuple with the frame counter and the frame itself.

Examples

julia> using MolSimToolkit, MolSimToolkit.Testing

julia> simulation = Simulation(
           Testing.namd_pdb, Testing.namd_traj; 
           first = 2, step = 2, last = 4
       );

julia> for frame in simulation 
           @show frame_index(simulation)
           # show x coordinate of first atom 
           @show positions(frame)[1].x
       end
frame_index(simulation) = 2
((positions(frame))[1]).x = 5.912472724914551
frame_index(simulation) = 4
((positions(frame))[1]).x = 7.346549034118652

julia> for (i, frame) in pairs(simulation)
           @show i, frame_index(simulation)
       end
(i, frame_index(simulation)) = (2, 2)
(i, frame_index(simulation)) = (4, 4)  

julia> for (i, frame) in enumerate(simulation)
           @show i, frame_index(simulation)
       end
(i, frame_index(simulation)) = (1, 2)
(i, frame_index(simulation)) = (2, 4)
source
Base.closeMethod
close(simulation::Simulation)

Closes the trajectory file.

source
Base.lengthMethod
length(simulation::Simulation)

Returns the number of frames to be iterated over in the trajectory file, considering the current frame range.

source
MolSimToolkit.first_frame!Method
first_frame!(simulation::Simulation)

Restarts the trajectory buffer, and places the current frame at the first frame in the trajectory.

Example

julia> using MolSimToolkit, MolSimToolkit.Testing

julia> simulation = Simulation(Testing.namd_pdb, Testing.namd_traj);

julia> first_frame!(simulation) 
Simulation 
    Atom type: Atom
    PDB file: structure.pdb
    Trajectory file: structure.dcd
    Total number of frames: 5
    Frame range: 1:1:5
    Number of frames in range: 5
    Current frame: 1
source
MolSimToolkit.frame_indexMethod
frame_index(simulation::Simulation)

Returns the index of the current frame in the trajectory. Returns nothing if no frame frame from the trajectory range has been read yet.

source
MolSimToolkit.get_frameMethod
get_frame(simulation::Simulation, iframe::Integer)

Returns the frame at the given index in the trajectory.

Example

julia> using MolSimToolkit, MolSimToolkit.Testing, PDBTools

julia> sim = Simulation(Testing.namd_pdb, Testing.namd_traj);

julia> frame4 = get_frame(sim, 4)
   Array{Atoms,1} with 20465 atoms with fields:
   index name resname chain   resnum  residue        x        y        z occup  beta model segname index_pdb
       1    N     ILE     P      211        1   -0.397   12.048   37.441  1.00  0.00     1    PROT         1
       2  HT1     ILE     P      211        1   -0.779   11.123   37.726  1.00  0.00     1    PROT         2
       3  HT2     ILE     P      211        1   -0.393   12.662   38.280  1.00  0.00     1    PROT         3
                                                       ⋮ 
   20463  SOD     SOD     S       13     4374  -11.686   23.749   19.935  1.00  0.00     1     SOD     20463
   20464  SOD     SOD     S       14     4375  -34.214   38.148   55.179  1.00  0.00     1     SOD     20464
   20465  SOD     SOD     S       15     4376    7.220  -52.702   66.223  1.00  0.00     1     SOD     20465

julia> writePDB(frame4, "frame4.pdb")
Note

The get_frame function will read the frames in the trajectory until the desired frame is reached. This can be slow for large trajectories. If the required frame is before the current frame of the simulation, the simulation will be restarted. The simulation object is returned positioned in the required frame.

source
MolSimToolkit.next_frame!Method
next_frame!(simulation::Simulation)

Reads the next frame in the trajectory file and returns it. Moves the current frame to the next one in the range to be considered (given by frame_range(simulation)).

source
MolSimToolkit.set_frame_range!Method
set_frame_range!(simulation::Simulation; first=1, last=nothing, step=1)

Resets the frame range to be iterated over. This function will restart the iteration of the simulation trajectory.

source

Positions

MolSimToolkit.positionsMethod
positions(frame::Chemfiles.Frame)

Return the positions of the atoms in a Chemfiles.Frame as a FramePositions object.

This is the default way to access the positions of the atoms in a simulation.

Example

julia> using MolSimToolkit, MolSimToolkit.Testing

julia> simulation = Simulation(Testing.namd_pdb, Testing.namd_traj);

julia> frame = current_frame(simulation);

julia> p = positions(frame);

julia> p[1].x 
5.912472724914551
source
MolSimToolkit.FramePositionsType
FramePositions{T,P<:Point3D{T},M<:AbstractArray{T}} <: AbstractVector{P}

Container for the positions of a set of atoms. The positions are stored in a matrix, where each column corresponds to the coordinates of an atom. The container is used such that using the coodinates from a Chemfiles.Frame is transparent to the user, and the coordinates can be accessed as p[i] where i is the index of the atom.

The coordinates of the atom can be accessed as p[i].x, p[i].y, and p[i].z.

A FramePositions object can be created with the positions function. The construction with FramePositions is not considered part of the public API.

Example

julia> using MolSimToolkit, MolSimToolkit.Testing

julia> simulation = Simulation(Testing.namd_pdb, Testing.namd_traj);

julia> frame = current_frame(simulation);

julia> p = positions(frame)
20465-element FramePositions{Float64, Point3D{Float64}, Chemfiles.ChemfilesArray}:
 [5.912472724914551, 10.768872261047363, 28.277008056640625]
 [5.040304183959961, 10.810898780822754, 27.71207046508789]
 ⋮
 [11.16289234161377, -37.30374526977539, 22.80788230895996]

julia> p[1]
3-element Point3D{Float64} with indices SOneTo(3):
  5.912472724914551
 10.768872261047363
 28.277008056640625

julia> p[1].x
5.912472724914551

julia> p[1].y
10.768872261047363

julia> p[1].z
28.277008056640625

It is also possible to take slices and views of the positions:

julia> p[1:2]
2-element FramePositions{Float64, Point3D{Float64}, Matrix{Float64}}:
 [5.912472724914551, 10.768872261047363, 28.277008056640625]
 [5.040304183959961, 10.810898780822754, 27.71207046508789]

julia> @view(coor[1:2])
 2-element FramePositions{Float64, Point3D{Float64}, SubArray{Float64, 2, Chemfiles.ChemfilesArray, Tuple{Base.Slice{Base.OneTo{Int64}}, UnitRange{Int64}}, false}}:
  [5.912472724914551, 10.768872261047363, 28.277008056640625]
  [5.040304183959961, 10.810898780822754, 27.71207046508789]
 
source

Unit cell

Wrap coordinates

MolSimToolkit.wrapFunction
wrap(x, xref, unit_cell_matrix::SMatrix{N,N,T}) where {N,T}
wrap(x, xref, sides::AbstractVector)

Wraps the coordinates of point x such that it is the minimum image relative to xref. The unit cell may be given a a static matrix of size (N,N) or as a vector of length N.

source
MolSimToolkit.wrap_to_firstFunction
wrap_to_first(x, unit_cell_matrix)

Wraps the coordinates of point x such that the returning coordinates are in the first unit cell with all-positive coordinates. The unit cell has to be a matrix of size (N,N).

Example

julia> using MolSimToolkit

julia> uc = [10 0 0; 0 10 0; 0 0 10]
UnitCell{Int64}
  10.000   0.000   0.000
   0.000  10.000   0.000
   0.000   0.000  10.000

julia> wrap_to_first(Point3D(15.0, 13.0, 2.0), uc)
3-element SVector{3, Float64} with indices SOneTo(3):
 5.0
 3.0000000000000004
 2.0
source