Ctrl+K

HF in VeloxChem

HF in VeloxChem#

import veloxchem as vlx

Define the molecule#

We first define the structure of a water molecule and choose a basis set.

mol_xyz = """3 

O       0.0000000000     0.1178336003     0.0000000000
H      -0.7595754146    -0.4713344012    -0.0000000000
H       0.7595754146    -0.4713344012     0.0000000000
"""

molecule = vlx.Molecule.read_xyz_string(mol_xyz)
basis = vlx.MolecularBasis.read(molecule, "cc-pVDZ", ostream=None)

molecule.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

print("Number of atoms:", molecule.number_of_atoms())
print("Number of electrons:", molecule.number_of_electrons())
print("Number of contracted basis functions:", basis.get_dimensions_of_basis(molecule))

Number of atoms: 3
Number of electrons: 10
Number of contracted basis functions: 24

SCF optimization#

Perform a self-consistent field (SCF) optimization to obtain the Hartree–Fock wave function and the associated ground-state energy.

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.ostream.mute()
scf_results = scf_drv.compute(molecule, basis)

SCF information#

The SCF driver object has a method named get_scf_energy() for retrieving the final energy.

print(f"Hartree–Fock energy: {scf_drv.get_scf_energy():14.10f} a.u.")

Hartree–Fock energy: -76.0265782198 a.u.

The return object from the compute() method is a Python dictionary containing several tensors:

  • C: molecular orbital coefficients as a NumPy array

  • E: orbital energies as a NumPy array

  • D: \(\alpha\)- and \(\beta\)-spin density matrices as a tuple of NumPy arrays

  • F: \(\alpha\)- and \(\beta\)-spin Fock matrices as a tuple of NumPy arrays

  • S: overlap integrals as a NumPy array

print("Dictionary keys:\n", scf_results.keys())
print("\nOrbital energies:\n", scf_results["E_alpha"])

Dictionary keys:
 dict_keys(['eri_thresh', 'qq_type', 'scf_type', 'scf_energy', 'restart', 'S', 'C_alpha', 'C_beta', 'E_alpha', 'E_beta', 'D_alpha', 'D_beta', 'F_alpha', 'F_beta', 'C', 'E', 'D', 'F'])

Orbital energies:
 [-20.55119961  -1.33473347  -0.69692651  -0.56605275  -0.49289827
   0.18486487   0.25566978   0.78604709   0.85103488   1.16387193
   1.20030818   1.25366153   1.44422793   1.4755769    1.67329727
   1.8679223    1.93111284   2.4430781    2.48048696   3.28268216
   3.334531     3.50569214   3.86065795   4.14355137]

Visualizing molecular orbitals#

The resulting molecular orbitals (MOs) can be visualized using the OrbitalViewer class.

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

In the eChem book, this is a static HTML page but if you download the notebook and open it in Jupyter lab you will be able to select and view any of the other molecular orbitals from the pull-down menu.

Contents