Periodic boundary conditions
MLatom currently supports calculations with periodic boundary conditions (PBC) only for ANI-type models (user-trained and universal models such as ANI-1ccx, ANI-2x but not their D4-corrected variants).
Warning
It is experimental option not extensively tested in simulations. No guarantees it works as you expect. We will be very thankful if you tell us about any problems and give suggestions by contacting us on Slack or by email.
Tutorial
Get started with examples on how to use it (notebook file):
import mlatom as ml
import numpy as np
Setting up the structure¶
We will use diamond to demonstrate in this tutorial.
Let's first define the structure of it.
a = 3.567 # the lattice constant of diamond
xyz = np.array(
[
[0, 0, 0],
[0, 2, 2],
[2, 0, 2],
[2, 2, 0],
[3, 3, 3],
[3, 1, 1],
[1, 3, 1],
[1, 1, 3]
]
) * a / 4 # the xyz coordinate for diamond crystal in the unicell
z = np.array([6] * 8) # eight carbons one cell...
mol = ml.molecule.from_numpy(xyz, z) # generate the molecule object
mol.pbc = True # equals to mol.pbc = [True, True, True]
mol.cell = a # equals to mol.cell = [a, a, a] or mol.cell = np.eye(3) *a
Check the structure with 3Dmol:
import py3Dmol
viewer = py3Dmol.view()
viewer.addModel(mol.get_xyz_string(), 'xyz')
viewer.setStyle({"stick":{"color": "white"}, "sphere": {"radius": 0.5, "color": "gray"}})
viewer.zoomTo()
viewer.show()
You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol
We can make it looks more diamonded with mlatom.data.molecule.proliferate():
import py3Dmol
viewer = py3Dmol.view()
viewer.addModel(mol.proliferate(XYZshifts=range(2)).get_xyz_string(), 'xyz')
viewer.setStyle({"stick":{"color": "white"}, "sphere": {"radius": 0.5, "color": "gray"}})
viewer.zoomTo()
viewer.show()
You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol
Choosing the model¶
model = ml.models.methods('ANI-1xnr')
Single-point calculations¶
After we defined our model, let's do single-point calculations
model.predict(molecule=mol,calculate_energy=True, calculate_energy_gradients=True)
# our energy should be -45.22730697839938 hartree
print(f'Energy: {mol.energy} hartree')
print('Energy gradients:')
print(mol.energy_gradients)
Energy: -45.22730697839938 hartree Energy gradients: [[-8.27383246e-08 -1.16474212e-07 -8.15539352e-08] [-7.74659838e-08 -2.06521236e-08 -5.23550625e-08] [-8.12822236e-08 -8.85738700e-08 -4.87361831e-08] [-8.77898856e-08 -1.75186869e-08 -6.22794687e-08] [ 1.23686505e-07 2.16007493e-07 2.48656761e-07] [ 1.18062871e-07 -6.88064574e-08 -1.25166935e-07] [ 3.73479452e-08 1.77298261e-07 -1.23865902e-07] [ 5.03121100e-08 -8.63303740e-08 2.39659130e-07]]
Geometry optimization¶
Now let's first add some noise to the geometry and then optimize it within the given cell.
Note: unfortunately, the cell parameters cannot be optimized yet.
mol_opt = ml.molecule.from_numpy(xyz + np.random.randn(*xyz.shape) * 0.1, z)
mol_opt.pbc = True
mol_opt.cell = a
ml.optimize_geometry(
model=model,
molecule=mol_opt,
program='ASE',
)
Step Time Energy fmax LBFGS: 0 02:46:37 -1225.742958 12.9199 LBFGS: 1 02:46:37 -1229.851896 4.3968 LBFGS: 2 02:46:37 -1230.601566 1.8222 LBFGS: 3 02:46:38 -1230.685981 0.5758 LBFGS: 4 02:46:38 -1230.697474 0.0887 LBFGS: 5 02:46:38 -1230.697700 0.0243 LBFGS: 6 02:46:38 -1230.697727 0.0024
<mlatom.simulations.optimize_geometry at 0x7f9b5deb47d0>
Compare geometries:
import py3Dmol
viewer = py3Dmol.view(width=960, height=320, viewergrid=(1, 3))
viewer.addModel(mol_opt.optimization_trajectory.to_database()[0].proliferate(XYZshifts=range(2)).get_xyz_string(), 'xyz', viewer=(0, 0))
viewer.addModel(mol_opt.proliferate(XYZshifts=range(2)).get_xyz_string(), viewer=(0, 1))
viewer.addModel(mol.proliferate(XYZshifts=range(2)).get_xyz_string(), viewer=(0, 2))
viewer.setStyle({"stick":{"color": "white"}, "sphere": {"radius": 0.5, "color": "gray"}})
viewer.addLabel("initial", {"alignment": "center", "screenOffset": {"x": 170, "y": -290}, "useScreen": True}, viewer=(0, 0))
viewer.addLabel("final", {"alignment": "center", "screenOffset": {"x": 170, "y": -290}, "useScreen": True}, viewer=(0, 1))
viewer.addLabel("reference", {"alignment": "center", "screenOffset": {"x": 170, "y": -290}, "useScreen": True}, viewer=(0, 2))
viewer.zoomTo()
viewer.show()
You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol
mol_md = ml.generate_initial_conditions(
molecule = mol_opt,
generation_method = 'random',
initial_temperature = 300,
random_seed = 0
)[0]
Start the dynamics:
dyn = ml.md(
model=model,
molecule=mol_md,
time_step=0.5,
ensemble='NVE',
maximum_propagation_time=32,
)
Visualize the trajectory:
import py3Dmol
viewer = py3Dmol.view()
viewer.addModelsAsFrames(dyn.molecular_trajectory.to_database().proliferate(XYZshifts=range(2)).get_xyz_string(), 'xyz',)
viewer.setStyle({"stick":{"color": "white"}, "sphere": {"radius": 0.5, "color": "gray"}})
viewer.zoomTo()
viewer.animate({'loop': "forward"})
viewer.show()
You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol
Limitations
Only ANI-type models:
universal ANI models ANI-1ccx, ANI-1xnr, ANI-2x
user-trained ANI-type machine learning potentials (only for making predictions)