System solvation

System solvation#

import veloxchem as vlx

Building systems#

VeloxChem can be used to solvate system in preparation for molecular dynamics simulations.

First, create a molecule object for the solute, in our example we choose deprotonated ibuprofen.

solute = vlx.Molecule.read_smiles("CC(C)CC1=CC=C(C=C1)C(C)C(=O)[O-]")

solute.show()

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

Second, generate a force field for the solute. By default RESP charges will be computed but, here, we instead make use of semi-empirical partial charges for faster execution.

ff_gen = vlx.MMForceFieldGenerator()

ff_gen.partial_charges = solute.get_partial_charges(solute.get_charge())
ff_gen.create_topology(solute)

Third, we use of the SolvationBuilder class to create a simulation box with a solute padding 1.0 nm (default) consisting of SPC/E water (default).

The solvator will run an NPT equilibration of 5 ps at 300K.

solvator = vlx.SolvationBuilder()

solvator.solvate(solute, equilibrate=True, neutralize=False)

solvator.system_molecule.show()

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

Mixed solvents#

Solvents can be mixed. Molecule objects are created for each type of solvent molecules.

water = vlx.Molecule.read_smiles("O")
propylene_glycol = vlx.Molecule.read_smiles("CC(CO)O")

Let us consider a solvation with 50/50 amounts of water and propylene glycol.

solvator.custom_solvate(
    solute,
    solvents=[water, propylene_glycol],
    proportion=(50, 50),
    box_size=(40, 40, 40),
)

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[8], line 1
----> 1 solvator.custom_solvate(
      2     solute,
      3     solvents=[water, propylene_glycol],
      4     proportion=(50, 50),
      5     box_size=(40, 40, 40),
      6 )

TypeError: SolvationBuilder.custom_solvate() got an unexpected keyword argument 'proportion'

solvator.system_molecule.show()

Free energy of solvation#

solvation = vlx.SolvationFepDriver()

# Settings are here chosen for a quick execution,
# it is recommended to use the default settings
solvation.num_steps = 5000
solvation.number_of_snapshots = 500

solvation_results = solvation.compute(solute, ff_gen)

print(
    f"Ibuprofen in water/propylene glycol: {solvation_results['free_energy']:.1f} kJ/mol"
)

Calculations based on files#

solvationfep = vlx.SolvationFepDriver()

# Using Gromacs type of input files
solvator.write_gromacs_files(solute_ff=ff_gen)
solvation_results = solvation.compute_solvation_from_gromacs_files(
    "system.gro", "system.top", "solute.gro", "solute.top"
)

# Using OpenMM type of input files
# With a mixed solvent there are multiple XML files in the list
solvator.write_openmm_files(solute_ff=ff_gen)
solvation_results = solvation.compute_solvation_from_omm_files(
    "system.pdb", "solute.pdb", "solute.xml", ["solvent_1.xml"]
)
Contents