Optimizing electronic structure#

Here we consider the optimization of the electronic structure using Hartree–Fock (HF) and Kohn-Sham density functional theory (KS-DFT), and how to plot the resulting molecular orbitals. Our example is methanol, using the xTB structure from the previous section.

import py3Dmol as p3d
import veloxchem as vlx

Specify structure#

The molecular geometry is written in xyz-coordinates, and visualized using py3Dmol.

methanol_xyz = """6

C        -0.362051846542         0.008519935441         0.024906559609
O         0.933468973663         0.509856441784        -0.193822841148
H        -0.483698601596        -0.368885697455         1.046180408743
H        -0.612499226929        -0.790059463155        -0.682044832355
H        -1.045997965440         0.841433017766        -0.127103979461
H         1.570778666846        -0.200864234382        -0.068115315388
"""

viewer = p3d.view(width=300, height=200)
viewer.addModel(methanol_xyz, 'xyz')
viewer.setViewStyle({"style": "outline", "width": 0.05})
viewer.setStyle({"stick":{},"sphere": {"scale":0.25}})
viewer.show()

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

SCF optimization#

At the HF/6-311G level of theory, the electronic structure is optimized by creating the molecule and basis objects, setting up an unrestricted SCF driver, and performing the SCF calculation:

molecule = vlx.Molecule.read_xyz_string(methanol_xyz)
basis = vlx.MolecularBasis.read(molecule, "6-311G")
scf_hf = vlx.ScfRestrictedDriver()
hf_results = scf_hf.compute(molecule, basis)
print(f'The total (HF) energy is {scf_hf.get_scf_energy(): .5f} Hartree')
The total (HF) energy is -115.01894 Hartree

Visualizing molecular orbitals#

The MOs can be visualized using the OrbitalViewer function:

viewer = vlx.OrbitalViewer()
viewer.plot(molecule, basis, scf_hf.mol_orbs)

Here you can only see the default MO selection (the HOMO), but if you run this in a notebook you can visualize any MO, change isovalues, and more.

KS-DFT optimization#

With the B3LYP exchange-correlation functional, the SCF is optimized using:

molecule = vlx.Molecule.read_xyz_string(methanol_xyz)
basis = vlx.MolecularBasis.read(molecule, "6-311G")
scf_b3lyp = vlx.ScfRestrictedDriver()
scf_b3lyp.xcfun = 'b3lyp'
b3lyp_results = scf_b3lyp.compute(molecule, basis)
print(f'The total (B3LYP) energy is {scf_b3lyp.get_scf_energy(): .5f} Hartree')
The total (B3LYP) energy is -115.71437 Hartree

And the MOs are visualized as:

viewer = vlx.OrbitalViewer()
viewer.plot(molecule, basis, scf_b3lyp.mol_orbs)