In [13]:
## Из-за того, что время от времени notebook перестает показывать выдачи, перезагружал регулярно
## Kernel
# Делаю тоннель
# plink -ssh -X iltarn@kodomo.cmm.msu.ru
# ipython notebook --no-browser
# plink -ssh -L 8864:localhost:8864 iltarn@kodomo.cmm.msu.ru
from xmlrpclib import ServerProxy
from IPython.display import Image
import os, sys
import scipy as sp
from scipy import constants
from scipy.constants import codata
import numpy as np
In [4]:
%%bash
export PATH=${PATH}:/home/preps/golovin/progs/bin
export MOPAC_LICENSE=/home/preps/golovin/progs/bin
export LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:/home/preps/golovin/progs/lib
In [7]:
%%bash
# Задание 1.
# Генерирую SMILES
cd
cd Term6
cd Practice2
echo 'c1cc2cc3ccc(cc4ccc(cc5ccc(cc1n2)[nH]5)n4)[nH]3' > porphyrin1.smi
In [9]:
%%bash
# Преобразую в .mol
obgen porphyrin.smi > porphyrin.mol
In [30]:
# Запускаю PyMol без графического окна
import __main__
__main__.pymol_argv = [ 'pymol', '-cp' ]
import pymol
pymol.finish_launching()
from pymol import cmd
In [31]:
from IPython.display import Image
In [14]:
# Устанавливаю настройки визуализации, подгружаю и визуализирую молекулу порфирина.
cmd.do('''
reinitialize
set ray_trace_mode, 0; red
set ray_opaque_background, off
set antialias, .5
set light_count, 8
set ambient, 0.5
set ray_trace_color, red
set cartoon_side_chain_helper, on
cd /home/students/y12/iltarn/Term6/Practice2
load porphyrin.mol
orient porphyrin
as sticks, porphyrin
zoom porphyrin, 4
label all, " %s" % ID
set label_color, black
ray
png pic1.png
''')
In [15]:
Image(filename='pic1.png')
Out[15]:
In [16]:
cmd.do('''
rotate x, 90
ray
png pic1_90.png
''')
In [17]:
Image(filename='pic1_90.png')
Out[17]:
In [18]:
# Исправляю молекулу и сохраняю в .pdb
cmd.do('''
remove id 39+40
save porphyrin.pdb
ray
png pic1_90_rem.png
''')
In [19]:
Image(filename='pic1_90_rem.png')
Out[19]:
In [345]:
%%bash
#С помощью openbabel форматирую координаты в mol формате во входной файл для Mopac.
# Используем метод параметризации PM6
babel -ipdb porphyrin.pdb -omop porphyrin_opt.mop -xk "PM6"
In [347]:
%%bash
# Скользкий шаг. Делал через Putty непосредственно - там все работает
#(правда, предупреждает о наличии новой версии и ждет нажатия Enter).
# На этом этапе происходит расчет молекулярных орбиталей
# для всех валентных электронов порфирина.
MOPAC2009.exe porphyrin_opt.mop
In [348]:
%%bash
# Переформатирование из .out в .pdb
babel -imopout porphyrin_opt.out -opdb porphyrin_opt.pdb
In [349]:
# Визуализация изменений.
cmd.do('''
delete all
load porphyrin_opt.pdb, porphyrin
orient porphyrin
as sticks, porphyrin
zoom porphyrin, 4
ray
png pic2.png
''')
In [350]:
# Наш порфирин стал плоским!
Image(filename='pic2.png')
Out[350]:
In [365]:
#Вид сбоку
cmd.do('''
rotate x, 90
ray
png pic2_90.png
''')
In [366]:
Image(filename='pic2_90.png')
Out[366]:
In [355]:
%%bash
# Сравним PM6 с AM1
babel -ipdb porphyrin.pdb -omop porphyrin_opt_am1.mop -xk "AM1"
MOPAC2009.exe porphyrin_opt_am1.mop
babel -imopout porphyrin_opt_am1.out -opdb porphyrin_opt_am1.pdb
In [369]:
cmd.do('''
delete all
load porphyrin_opt_am1.pdb, porphyrin_am1
load porphyrin_opt.pdb, porphyrin_pm6
color red, porphyrin_am1
color blue, porphyrin_pm6
orient porphyrin
as sticks, porphyrin_pm6
as sticks, porphyrin_am1
zoom porphyrin, 4
rotate x, 30
ray
png pic3.png
''')
In [370]:
Image(filename='pic3.png')
# Видим, что на структурах сильнее расходятся водороды.
Out[370]:
In [379]:
%%bash
# Задание 2. Рассчет возбужденных состояний порфирина
# Копирую файл .mop
*porphyrin_opt.mop* porphyrin_opt_spectr.mop
echo "" >> porphyrin_opt_spectr.mop
echo "cis c.i.=4 meci oldgeo" >> porphyrin_opt_spectr.mop
echo "some description" >> porphyrin_opt_spectr.mop
MOPAC2009.exe porphyrin_opt_spectr.mop
In [22]:
#посмотрим на файл
%load /home/students/y12/iltarn/Term6/Practice2/porphyrin_opt_spectr.out
In [11]:
# После просмотра файла выведем отдельно энергии электронных переходов ( в самом конце).
f=open('/home/students/y12/iltarn/Term6/Practice2/porphyrin_opt_spectr.out','r')
lines=f.readlines()
f.close()
for l in lines[-19:-6]:
print l.split()
In [23]:
# Рассчитаем спектр поглощения
for i in lines[-15:-7]:
l=i.split()
k=float(l[2])
print l[0], l[4], (constants.c*constants.h*(10**9)/codata.value('electron volt-joule relationship'))/k, 'nm'
f.close()
In []:
# Видим, что порфирин поглащает в фиолетовой области спектра , а также в видимой.
In []:
# Задание 3. Определение геометрии молекулы O=C1C=CC(=O)C=C1.
In [28]:
%%bash
cd
cd Term6/Practice2
# Сгенерируем файл .smi, затем построим по нему 3D.
echo 'O=C1C=CC(=O)C=C1 benzoquinone' > benzoquinone.smi
obgen benzoquinone.smi > benzoquinone.mol
In [47]:
# Устанавливаю настройки визуализации, подгружаю и визуализирую молекулу бензохинона.
cmd.do('''
reinitialize
set ray_trace_mode, 0; red
set ray_opaque_background, off
set antialias, .5
set light_count, 8
set ambient, 0.5
set ray_trace_color, red
set cartoon_side_chain_helper, on
cd /home/students/y12/iltarn/Term6/Practice2
load benzoquinone.mol, benzoquinone
orient benzoquinone
as sticks, benzoquinone
zoom benzoquinone, 1
ray
png pic4.png
''')
In [48]:
Image(filename='pic4.png')
Out[48]:
In [60]:
#Вид сбоку
cmd.do('''
rotate x, -70
ray
png pic4_70.png
save benzoquinone.pdb
''')
In [58]:
Image(filename='pic4_70.png')
Out[58]:
In [51]:
# Видно, что молекула плоская.
In [62]:
%%bash
# Запустим MOPAC
babel -ipdb benzoquinone.pdb -omop benzoquinone.mop -xk "PM6"
MOPAC2009.exe benzoquinone.mop
babel -imopout benzoquinone.out -opdb benzoquinone_mopac.pdb
In [78]:
# Сравним молекулы:
cmd.do('''
delete all
load benzoquinone_mopac.pdb, benzoquinone_mopac
load benzoquinone.pdb, benzoquinone
color red, benzoquinone
color blue, benzoquinone_mopac
set valence, on
rotate x, -60
ray
png pic5.png
''')
In [79]:
Image(filename='pic5.png')
Out[79]:
In []:
# После MOPAC молекула осталась плоской, однако некоторые небольшие различия все же присутствуют.
# Посмотрим, что будет в случае дианиона...
In [80]:
with open("benzoquinone.mop") as f,open("benzoquinone_charge.mop",'w') as o:
data=f.read()
data=data.replace("PM6\benzoquinone.pdb","PM6 CHARGE=-2\benzoquinone.pdb")
o.write(data)
with open("benzoquinone_charge.mop") as f,open("benzoquinone_charge2.mop",'w') as o:
data=f.read()
data=data.replace("O","O(-)")
o.write(data)
In [81]:
%%bash
MOPAC2009.exe benzoquinone_charge2.mop
babel -imopout benzoquinone_charge2.out -opdb benzoquinone_charge2.pdb
In [82]:
cmd.do('''
delete all
load benzoquinone_charge2.pdb, anion
load benzoquinone.pdb, benzoquinone
color red, benzoquinone
color blue, anion
set valence, on
rotate x, -60
ray
png pic6.png
''')
In [83]:
Image(filename='pic6.png')
Out[83]:
In []:
# В анионе видно, что атомы кислорода сместились в стороны друг от друга (синий).
#В целом, судя по синему контуру снаружи кольца, кольцо подрастянулось.
In [92]:
# Задание 4. Тиминовые димеры.
cmd.do('''
delete all
load http://kodomo.fbb.msu.ru/FBB/year_08/term6/td.pdb
as sticks, td
set valence, on
ray
png pic7.png
''')
Image(filename='pic7.png')
Out[92]:
In [95]:
%%bash
# Сначала оптимизируем структуру с PM6 с зарядом 0.
wget "http://kodomo.fbb.msu.ru/FBB/year_08/term6/td.pdb" -O td.pdb
babel -ipdb td.pdb -omop td.mop -xk "PM6"
MOPAC2009.exe td.mop
babel -imopout td.out -opdb td_opt.pdb
In [110]:
%%bash
# Теперь с зарядом +2
babel -ipdb td_opt_2.pdb -omop td_opt_2.mop -xk "PM6"
sed -i 's/PM6/PM6 CHARGE=+2/' td_opt_2.mop
In [113]:
%%bash
MOPAC2009.exe td_opt_2.mop
babel -imopout td_opt_2.out -opdb td_opt_2.pdb
In [120]:
cmd.do('''
delete all
load td_opt_2.pdb
as sticks, td_opt_2
set valence, on
ray
png pic8.png
''')
Image(filename='pic8.png')
Out[120]:
In [123]:
%%bash
# Возвращаем 0 заряд:
babel -ipdb td_opt_2.pdb -omop td_opt_0.mop -xk "PM6"
sed -i 's/PM6 CHARGE=+2/PM6/' td_opt_0.mop
MOPAC2009.exe td_opt_0.mop
babel -imopout td_opt_0.out -opdb td_opt_0.pdb
In [127]:
cmd.do('''
delete all
load td_opt_0.pdb
as sticks, td_opt_0
set valence, on
ray
png pic9.png
''')
In [128]:
Image(filename='pic9.png')
Out[128]:
In []:
# Тиминовый димер не приобрел той же структуры, что и при исходно 0 заряде.
In [24]:
# Сравним энергии... Строки подбирались после просмотра файлов.
f=open('/home/students/y12/iltarn/Term6/Practice2/td.out','r')
lines=f.readlines()
f.close()
print lines[205]
f=open('/home/students/y12/iltarn/Term6/Practice2/td_opt_0.out','r')
lines=f.readlines()
f.close()
print lines[364]
f=open('/home/students/y12/iltarn/Term6/Practice2/td_opt_2.out','r')
lines=f.readlines()
f.close()
print lines[277]
In [25]:
# Как видно, состояния тиминового димера (-3273.58217 EV) и двух несвязанных тиминовых нуклеотидов
# (-3273.69661 EV) примерно одинаковы по энергии, однако второе состояние более предпочтительно
# ввиду большей (по модулю) энегрии. Поэтому возбужденное состояние и свалилось в отдельные тимины.