Solvate#

Adding solvent molecules and ions to solvate a molecular system.

The function molsysmt.build.solvate() embeds a solute within a solvent box (such as explicit TIP3P water), adjusts periodic boundary dimensions, and adds neutralizing and physiological ions.

Added in version 1.0.0.

API documentation

Follow this link for a detailed description of the input arguments, raised errors, and returned objects of this function: molsysmt.build.solvate().

Basic usage#

Let’s show how to solvate a protein structure loaded from PDB ID 1VII:

import molsysmt as msm
import numpy as np

molsys = msm.convert('pdb_id:1vii', to_form='molsysmt.MolSys')
molsys = msm.basic.remove(molsys, selection='molecule_type=="water"')

We check whether the initial system is solvated using molsysmt.build.is_solvated():

msm.build.is_solvated(molsys)

False

Now we solvate the system in a cubic box with 14 Å clearance using molsysmt.build.solvate(). Solvation requires exactly one structure because independently solvated frames can produce different solvent topologies. For an ensemble or trajectory, first select one representative structure with msm.convert(..., structure_indices=index):

molsys_cub = msm.build.solvate(molsys, box_shape='cubic', clearance='14.0 angstroms')

msm.build.is_solvated(molsys_cub)

True

We inspect the box vectors and angles of the solvated system using molsysmt.basic.get():

box, box_angles, box_shape = msm.get(molsys_cub, element='system', box=True, box_angles=True, box_shape=True)
print('Box:', box)
print('Angles:', box_angles)
print('Shape:', box_shape)

Box: [[[4.930948498633193 0.0 0.0]  [0.0 4.930948498633193 0.0]  [0.0 0.0 4.930948498633193]]] nanometer
Angles: [[1.570796 1.570796 1.570796]] radian
Shape: cubic

We visualize the solvated system with molsysmt.basic.view().

solvate returns the system in MolSysMT’s periodic box convention: the box vectors describe shape and size, and the cell they span starts at the origin. Molecules are kept whole, so the few that cross a face have some atoms sitting just outside the drawn cell. That is the correct picture of a periodic system, not a defect: the cell tiles space, and a molecule crossing a wall reappears on the opposite side.

Wrapping every atom into the cell would remove those overhangs and split those molecules instead, drawing bonds that stretch across the whole box. If you want that — because a simulation engine is going to re-image anyway — ask for it explicitly with molsysmt.pbc.wrap_to_pbc() and compact=False.

msm.view(molsys_cub)

Adding physiological ions#

Passing ionic_strength (e.g. '150.0 millimolar') automatically introduces neutralizing counterions and physiological salt concentration (NaCl):

molsys_ions = msm.build.solvate(molsys, box_shape='cubic', clearance='14.0 angstroms', ionic_strength='150.0 millimolar')
n_waters = msm.get(molsys_ions, element='system', n_waters=True)
n_ions = msm.get(molsys_ions, element='system', n_ions=True)
print(f'Solvated system contains {n_waters} water molecules and {n_ions} ions.')

Solvated system contains 3545 water molecules and 22 ions.

Box shape#

MolSysMT supports multiple periodic box geometries:

'cubic': Standard rectangular/cubic unit cell.
'truncated octahedral': Truncated octahedron unit cell, which reduces the number of solvent molecules required to maintain a given clearance.
'rhombic dodecahedral': Rhombic dodecahedron unit cell, another space-filling lattice optimal for isotropic globular solutes.

Let’s compare the water requirements across these different geometries for the same clearance:

molsys_oct = msm.build.solvate(molsys, box_shape='truncated octahedral', clearance='14.0 angstroms')
molsys_rhomb = msm.build.solvate(molsys, box_shape='rhombic dodecahedral', clearance='14.0 angstroms')

n_waters_cub = msm.get(molsys_cub, element='system', n_waters=True)
n_waters_oct = msm.get(molsys_oct, element='system', n_waters=True)
n_waters_rhomb = msm.get(molsys_rhomb, element='system', n_waters=True)

print(f'Cubic box: {n_waters_cub} waters (100% reference)')
print(f'Truncated octahedral box: {n_waters_oct} waters (saved {100.0 * (n_waters_cub - n_waters_oct) / n_waters_cub:.1f}%)')
print(f'Rhombic dodecahedral box: {n_waters_rhomb} waters (saved {100.0 * (n_waters_cub - n_waters_rhomb) / n_waters_cub:.1f}%)')

Cubic box: 3565 waters (100% reference)
Truncated octahedral box: 1454 waters (saved 59.2%)
Rhombic dodecahedral box: 1299 waters (saved 63.6%)

We visualize the truncated octahedral solvated system. The same convention applies: whole molecules, cell starting at the origin, and a thin shell of atoms beyond the faces.

msm.view(molsys_oct)