Kodomo

Пользователь

Учебная страница курса биоинформатики,
год поступления 2013

Вычисление атомных орбиталей водорода

Цель занятия: опираясь на уравнение построить электронную плотность в одно-электронном атоме и сравнить с известными программами

Работа разделена на две части: расчёт плотности в Ipython Notebook с сохранением в формате CUBE и визуализация в Pymol.

   1 import numpy
   2 import scipy.special
   3 import scipy.misc

   1 import npy2cube

   1 def w(n,l,m,d):
   2     
   3     x,y,z = numpy.mgrid[-d:d:30j,-d:d:30j,-d:d:30j]
   4     
   5     r = lambda x,y,z: numpy.sqrt(x**2+y**2+z**2)
   6     theta = lambda x,y,z: numpy.arccos(z/r(x,y,z))
   7     phi = lambda x,y,z: numpy.arctan(y/x)
   8     
   9     a0 = 1.
  10     
  11     R = lambda r,n,l: (2*r/n/a0)**l * numpy.exp(-r/n/a0) * scipy.special.genlaguerre(n-l-1,2*l+1)(2*r/n/a0)
  12     WF = lambda r,theta,phi,n,l,m: R(r,n,l) * scipy.special.sph_harm(m,l,phi,theta)
  13     absWF = lambda r,theta,phi,n,l,m: numpy.absolute(WF(r,theta,phi,n,l,m))**2
  14       
  15     return WF(r(x,y,z),theta(x,y,z),phi(x,y,z),n,l,m)

   1 # определите шаг grid при заданном диапозоне от -d до d
   2 d=.....
   3 step= .....
   4 
   5 # Зададим цикл по перебору квантовых чисел
   6 
   7 for n in range(0,4):
   8     for l in range(0,n):
   9         for m in range(0,l+1,1):
  10             grid= .... 
  11             name='%s-%s-%s' % (n,l,m)
  12             # для сохранения нужно задать координаты старта grid и шаг по каждому направлению
  13             npy2cube.npy2cube(grid,(-d,....),(step,.....),name+'.cube')

   1 ### Откуда эти цифры? 
   2 pml='''
   3 ### cut below here and paste into script ###
   4 cmd.volume_ramp_new('ramp007', [\
   5       0.005, 0.00, 0.00, 1.00, 0.00, \
   6       0.01,  0.00, 1.00, 1.00, 0.20, \
   7       0.015, 0.00, 0.00, 1.00, 0.00, \
   8     ])
   9 ### cut above here and paste into script ###
  10 '''
  11 
  12 for n in range(0,4):
  13     for l in range(0,n):
  14         for m in range(0,l+1,1):
  15             name='%s-%s-%s' % (n,l,m)
  16             # Загрузить cube файл
  17             # Отобразить электронную плотность (hint: volume)
  18             # Покрасить ее разумным образом 
  19 
  20 # запишите переменную в файл
  21 v=open(.....

   1 from IPython.display import display,Image
   2 display(Image=(file))

   1 export PATH=${PATH}:/tmp
   2 orca h.inp