Generating numpy-compatible electron density from Molden / QC calculation

Question

I'd like to use a variety of Python tools (e.g., matplotlib and Mayavi) to analyze orbitals plot electron density contours, etc. Most of these expect to work on numpy matrices.

My ideal workflow would be something like this:

import packageX

molecule = packageX.read('file.molden')
# specify grid options
density = packageX.calculateDensity(molecule)

Are there Python packages that make it easy to calculate electron density, spin density, and orbitals on a NumPy grid?

Answer 1

ORBKIT

This uses the ORBKIT package.

Gunter Hermann, Vincent Pohl, Jean Christophe Tremblay, Beate Paulus, Hans-Christian Hege, and Axel Schild, "ORBKIT: A Modular Python Toolbox for Cross-Platform Postprocessing of Quantum Chemical Wavefunction Data", J. Comput. Chem. 2016, 37, 1511-1520.

import orbkit as ok

ok.options.filename = '/path/to/86.molden' # also directly reads file.molden.gz if needed
ok.options.itype = 'molden'
ok.options.no_output = True # we'll directly get a Numpy array without generating a file
ok.options.adjust_grid = [5, 0.1] # 5 bohr padding on all sides, 0.1 bohr resolution
ok.options.numproc = 4 # or whatever you want - uses Python subprocess

data = ok.run_orbkit() # prints useful log bits in Jupyter

data.shape # returns (186, 179, 282) for my example file

If you want to do this in a loop, you'll need to set the ok.options every time:

ok.init()
# Set some options (again)
ok.options.filename = '/path/to/new/file.molden'
ok.options.itype = 'molden'
ok.options.no_output = True
ok.options.adjust_grid = [5, 0.1] # 5 bohr padding on all sides, 0.1 bohr resolution
ok.options.numproc = 4 # or whatever you want - uses Python subprocess
# run orbkit again
data_2 = ok.run_orbkit()

Orbkit easily handles the following formats:

• Molden
• wfn/wfx
• fchk
• Interface to cclib parser (for other formats)

Beyond electron density it's pretty easy to generate:

• orbitals
• spin density
• derivatives (analytical first and second order)
• Laplacian

More can be found on the ORBKIT website

I had to patch ORBKIT to read the Molden files I had (https://github.com/orbkit/orbkit/pull/39) but you can install via pip install orbkit or python setup.py install once you download.