Задание 7. Изучение работы методов контроля температуры в молекулярной динамике¶
Подготовительные работы¶
Импортируем нужные модули.
In [1]:
from subprocess import run
Скачаем в рабочую директорию все необходимые файлы.
In [2]:
%%bash
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/etane.gro -O et.gro
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/be.mdp
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/vr.mdp
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/nh.mdp
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/an.mdp
wget http://kodomo.fbb.msu.ru/FBB/year_08/term6/sd.mdp
Теперь допишем заготовку файла с топологией этана. Во-первых, дописали все необходимые связи, валентные и торсионные углы. Во-вторых, поменяли типы атомов: CX на opls_135, а HX на opls_140. В-третьих, определили, что тип потенциала торсионных углов равен 3. Итоговый файл выглядит так:
In [3]:
%%bash
cat et.top
И создадим файл для подсчёта длины C-C связи.
In [4]:
%%bash
echo -e "[ b ]\n1 2" > b.ndx
In [5]:
%%bash
cat b.ndx
Расчёт работы термостатов¶
Посчитаем всё, что от нас хотят, в цикле.
In [2]:
thermo_types = ['be', 'vr', 'nh', 'an', 'sd']
In [8]:
for thermo in thermo_types:
# configs
run(f'grompp -f {thermo}.mdp -c et.gro -p et.top -o et_{thermo}.tpr', shell=True)
print(thermo + '\tconfig')
# mdrun
run(f'mdrun -deffnm et_{thermo} -v -nt 1', shell=True)
print(thermo + '\tmdrun')
# pdb
run(f'echo 0 | trjconv -f et_{thermo}.trr -s et_{thermo}.tpr -o et_{thermo}.pdb', shell=True)
print(thermo + '\tpdb')
# energies
run(f'echo 10 11 | g_energy -f et_{thermo}.edr -o et_{thermo}_en.xvg -xvg none', shell=True)
print(thermo + '\tenergies')
# C-C bond lengths
run(f'g_bond -f et_{thermo}.trr -s et_{thermo}.tpr -o bond_{thermo}.xvg -n b.ndx -xvg none', shell=True)
print(thermo + '\tbond lengths')
Визуальный анализ¶
Загрузим pymol. Он не заработал на kodomo, поэтому переехал к себе на ноутбук.
In [3]:
import __main__
__main__.pymol_argv = [ 'pymol', '-x' ]
import pymol
pymol.finish_launching()
from pymol import cmd
In [6]:
import os
import imageio
import time
Всего было 251 шагов.
In [11]:
for thermo in thermo_types:
# aesthetics
cmd.reinitialize()
cmd.set('sphere_scale', 0.2)
cmd.set('stick_radius', 0.1)
cmd.set('ray_trace_mode', 3)
cmd.set('antialias', 2)
cmd.set('sphere_scale', 0.2)
cmd.set('stick_radius', 0.1)
# load data
cmd.load(f'et_{thermo}.pdb')
cmd.show('spheres')
cmd.orient()
# save images
os.mkdir(thermo)
cmd.mpng(os.path.join(thermo, ''), mode=2, width=640, height=480)
while not os.path.exists(os.path.join(thermo, '0251.png')):
time.sleep(0.1)
# create gif
with imageio.get_writer(f'{thermo}.gif', mode='I', fps=30) as writer:
for i in range(0, 251):
image = imageio.imread(os.path.join(thermo, f'{i + 1:04}.png'))
writer.append_data(image)
Теперь посмотрим на гифки.
In [12]:
from IPython.display import display, Image, HTML
In [14]:
thermo_names = ['Berendsen', 'V-rescale', 'Nose-Hoover', 'Andersen', 'Stochastic']
In [15]:
for thermo, name in zip(thermo_types, thermo_names):
display(HTML(name))
display(Image(f'{thermo}.gif', format='png'))
