Psi4 및 Python으로 프론티어 궤도 시각화

Jan 15 2023
Frontier Orbital은 HOMO(Highest Occupied Molecular Orbital) 및 LUMO(Lowest Occupied Orbital)를 의미합니다. 대부분의 유기 화학자들은 주로 이러한 오비탈에만 관심이 있는데, 이는 높은 에너지를 가진 전자가 다음 오비탈로 점프할 수 있는 상호 작용 지점이기 때문입니다.

Frontier Orbital은 HOMO(Highest Occupied Molecular Orbital) 및 LUMO(Lowest Occupied Orbital)를 의미합니다. 대부분의 유기 화학자들은 주로 이러한 오비탈에만 관심이 있는데, 이는 높은 에너지를 가진 전자가 다음 오비탈로 점프할 수 있는 상호 작용 지점이기 때문입니다. 그것들은 서로 가장 가능성이 높은 분자의 지배적인 상호 작용입니다. HOMO/LUMO 다이어그램의 예가 그림 1 에 나와 있습니다.

그림 1: HOMO/LUMO 분자 오비탈의 예

그러나 우리는 전자가 파동이고 이러한 오비탈에서 생성되는 전자 구름이 있다는 것을 알고 있습니다. 이 구름을 시각화하면 전자가 어떻게 행동하는지에 대한 아이디어와 몇 가지를 추론할 수 있는 기억에 남는 모양을 얻을 수 있습니다. 이를 위해 이전 블로그 게시물 에서 볼 수 있는 psi4 메탄올 z-매트릭스로 돌아가 보겠습니다 . 코드를 바로 살펴보겠습니다.

import psi4

psi4.set_memory('1000mb')
psi4.core.set_num_threads(1)

zmatrix = '''\
    H11
    O11  H11  0.9600
    C11  O11  1.4000  H11  108.0000
    H12  C11  1.1000  O11  112.0000  H11  -60.0000
    H13  C11  1.1000  O11  112.0000  H11   60.0000
    H14  C11  1.1000  O11  111.8699  H11 -180.0000
'''

universe = psi4.geometry(zmatrix)
universe.update_geometry()
universe.print_out()

energy, wave_function = psi4.optimize(
  'hf/6-31g*',
   return_wfn =True,
  molecule=universe
)

그래서 우리는 어떤 오비탈이 HOMO/LUMO 쌍인지 결정해야 하기 때문에 두 가지, 즉 알파 전자의 수와 분자 오비탈의 수를 원합니다. Pauli Exclusion Principle을 준수하기 때문에 2를 사용합니다.

number_of_occupied = 2
number_of_virutal = 2

alpha_electrons = wave_function.nalpha()
molecule = wave_function.nmo()

min_orb = max(1, alpha_electrons + 1 - number_of_occupied)
max_orb = min(molecule, alpha_electrons + number_of_virutal)
orbitals = [k for k in range(min_orb, max_orb + 1)]

cubeprop_tasks = []
  cubeprop_tasks.append('DENSITY')
  cubeprop_tasks.append('ORBITALS')

  psi4.set_options(
    {
      'scf_type': 'df',
      'g_convergence': 'gau_tight',
      'freeze_core': 'true',
      'CUBEPROP_TASKS': cubeprop_tasks,
      'CUBEPROP_FILEPATH': './',
      'CUBEPROP_ORBITALS': orbitals,
    }
  )

  psi4.cubeprop(wave_function)

Da.cube
Psi_a_9_2-A.cube

pip install moly
pip install opt_einsum

import moly

fig = moly.Figure()
fig.add_cube(
    'Psi_a_9_2-A.cube',
    iso=0.03,
    colorscale="rdbu",
    opacity=0.1
)

fig.show()