1) модуль pybel
import pybel
mol=pybel.readstring('smi','c1cc2cc3ccc(cc4ccc(cc5ccc(cc1n2)[nH]5)n4)[nH]3')
mol.addh()
mol.make3D(steps=10000)
mol.write(format='pdb',filename='pr5_1.pdb',overwrite=True)
Неправильная структура: не плоская и водороды в центре.
2). MOPAC способ PM6/АМ1
mop=mol.write(format='mopin',opt={'k':'PM6 CHARGE=%d' % mol.charge})
mop=mol.write(format='mopin',filename='my.mop',opt={'k':'PM6 CHARGE=%d' % mol.charge},overwrite=True)
В my.mop заряд молекулы равен нулю, поэтому перейдем к способу АМ1
mop_am1=mol.write(format='mopin',opt={'k':'AM1 CHARGE=%d' % mol.charge})
mop_am1=mol.write(format='mopin',opt={'k':'AM1 CHARGE=%d' % mol.charge},overwrite=True)
print mop_am1
%%bash
export MOPAC_LICENSE='/home/preps/golovin/progs/mopac/'
/home/preps/golovin/progs/mopac/MOPAC2016.exe my.mop
opt=pybel.readfile('mopout','my.out')
for i in opt:
print i
i.write(format='pdb',filename='my.pdb')
import subprocess
def mopac_pr5(smiles, method, prefix):
#стоило как-то иначе назвать функцию, но не хочется все переписывать
mol=pybel.readstring('smi',smiles)
mol.addh()
mol.make3D(steps=10000)
mop=mol.write(format='mopin',opt={'k': method+' CHARGE=%d' % mol.charge})
mop=mol.write(format='mopin',filename = prefix+'.mop',opt={'k': method+' CHARGE=%d' % mol.charge},overwrite=True)
mol.write(format='pdb',filename=prefix+'_before.pdb',overwrite=True)
cmd1 = " export MOPAC_LICENSE='/home/preps/golovin/progs/mopac/' "
cmd2 = '/home/preps/golovin/progs/mopac/MOPAC2016.exe '+prefix+'.mop'
subprocess.call(cmd1, shell=True)
subprocess.call(cmd2, shell=True)
opt=pybel.readfile('mopout',prefix+'.out')
for i in opt:
print i
i.write(format='pdb',filename=prefix+'.pdb',overwrite=True)
mopac_pr5('CC(=O)C','PM6', 'acetone')
mopac_pr5('C1=CC2=CC3=CC=C(N3)C=C4C=CC(=N4)C=C5C=CC(=N5)C=C1N2','PM6', 'por_pm6')
mopac_pr5('C1=CC2=CC3=CC=C(N3)C=C4C=CC(=N4)C=C5C=CC(=N5)C=C1N2','AM1', 'por_am1')
Методом AM1 не получился плоский порфирин, а PM6 получился. На рисунке - наложение структур, беленький - не плоский.
%%bash
cp my.mop my_spectr.mop
echo "" >> my_spectr.mop
echo "cis c.i.=4 meci oldgeo" >> my_spectr.mop
echo "some description" >> my_spectr.mop
%%bash
export MOPAC_LICENSE='/home/preps/golovin/progs/mopac/'
/home/preps/golovin/progs/mopac/MOPAC2016.exe my_spectr.mop
opt=pybel.readfile('mopout','my_spectr.out')
for i in opt:
print i
i.write(format='pdb',filename='my_spectr.pdb',overwrite=True)
from IPython.core.display import HTML
f = open('my_spectr.out', 'r')
lines = f.readlines()[-21:-13]
f.close()
energies = [float(line.split()[1]) for line in lines]
print energies, 'eV'
wavelengths = [1239.84193/e for e in energies]
frequencies = [e/4.135667516 for e in energies]
print 'frequncy,THz\t wavelength,nm\tenergy,eV'
for i in range(len(energies)):
print '%.4f \t %.4f \t %.4f' % (frequencies[i]*1000, wavelengths[i], energies[i])
Тимидиновые димеры
import pybel
mol=pybel.readfile("pdb", "td.pdb").next()
mop=mol.write(format='mop',filename='td.mop',opt={'k':'PM6 CHARGE=%d' % mol.charge},overwrite=True)
%%bash
export MOPAC_LICENSE='/home/preps/golovin/progs/mopac/'
/home/preps/golovin/progs/mopac/MOPAC2016.exe td.mop
opt=pybel.readfile('mopout','td.out')
for i in opt:
print i
i.write(format='pdb',filename='td_opt.pdb',overwrite=True)
Исходная и оптимизированная (розовенькая) структуры.
%%bash
cp td.mop td2.mop
#change line PM6 CHARGE=0 to PM6 CHARGE=+2
! /home/preps/golovin/progs/mopac/MOPAC2016.exe td2.mop
opt=pybel.readfile('mopout','td2.out')
for i in opt:
print i
i.write(format='pdb',filename='td2_opt.pdb',overwrite=True)
Два тимина, соединены все равно двумя связями.
mol=pybel.readfile("pdb", "td2_opt.pdb").next()
mop=mol.write(format='mop',filename='td2_opt.mop',opt={'k':'PM6 CHARGE=%d' % mol.charge},overwrite=True)
! /home/preps/golovin/progs/mopac/MOPAC2016.exe td2_opt.mop
opt=pybel.readfile('mopout','td2_opt.out')
for i in opt:
print i
i.write(format='pdb',filename='td2_last.pdb',overwrite=True)
Оптимизировали их:
%%bash
echo dimer
grep 'TOTAL ENERGY' td.out
echo semi-dimer
grep 'TOTAL ENERGY' td2.out
echo two thymins
grep 'TOTAL ENERGY' td2_opt.out
Ab initio вычисления: Gamess US
mopac_pr5('C1=CC=C2C=CC=CC2=C1','PM6','na')
mopac_pr5('C1=CC=C2C=CC=C2C=C1','PM6','az')
Азален до и после оптимизации (где какой и так понятно):
Нафтален до (светло-розовый) и после оптимизации:
mopac_pr5('c1ccc2ccccc2c1','PM6','na_ar')
mopac_pr5('c1ccc2cccc2cc1','PM6','az_ar')
Ароматический азулен до и после оптимизации: Ароматический нафтален до и после оптимизации (фиолетовый):
opt=pybel.readfile('mopout','na_ar.out')
for i in opt:
print i
i.write(format='orcainp',filename='orca.inp',overwrite=True)
%%bash
/srv/databases/orca/orca orca.inp
Умываю руки.