Barely Functional Theories

Musings on science and game design by James Furness.



Vibrational Analysis



A python script to solve the nuclear Schrödinger equation in the Born-Oppenheimer approximation for a diatomic molecule. The script comes in two files, the actual integrator Solver.py, the simple GUI wrapper QuantumWobbler.py and an example script for plugging the output of the solver into matplotlib. All three are provided here along with sample input, packaged as a zip archive shared here under the MIT licence.

Download QuantumWobbler.zip


The first 10 vibrational energy levels and wave functions for molecular hydrogen.

The first 10 vibrational energy levels and wave functions for molecular hydrogen.


About

This solver was written as the start of a side project to investigate the vibrational spectra of diatomic molecules in strong magnetic fields bound by the perpendicular paramagnetic bonding mechanism[1,2,3]. This project was sidelined before any magnetic fields could be incorporated, though the field free solver was finished. At its core, the script is simply an implementation of the Cooley-Numerov integrator [4] that solves the nuclear Schrödinger for a given guess energy. The output was polished up for use in a talk and the solver is uploaded here in case anyone would like to use it.

A GUI was added for use in undergraduate labs, but was not necessary in the end. As such it remains a little underdeveloped and buggy in places, though the core solver is solid. Whilst the GUI runs acceptably, it is quite unstable and not recommended for routine use. If anyone would like to use it for teaching please let me know and I will clean it up, otherwise it doesn’t warrant the effort to polish something that will not get used.

A much more refined solver called LEVEL has been developed by R. J. Le Roy, and has a far greater feature set. It is an old program however, and can be a little un-intuitive to use. If your needs are simple energy levels or a plot of the nuclear wave functions however, then my solver might fit your needs with less hassle.

The script is written in Python 2.7 and the core solver requires only the matplotlib and scipy modules to run. The GUI also requires the enthought.traits, chaco and pyface modules to run. The simplest way to acquire and install these modules is to use the package manager in a Python distribution such as Enthought Canopy, setup should be straight forwards.


General Operation

An input file containing information about the system is required for both GUI operation and when the solver is used directly.
It is essential for the solver’s operation that the sampled points on the electronic energy curve have a regular separation, e.g. 0.02,0.04,0.06… The solver is not sophisticated enough to interpolate between the data points provided.

Guess levels are optional and not required. If manual guesses are missing a fitted Morse potential will be used to generate automatic guess levels. Manually setting initial guesses can be useful when the automatic guess converges to the incorrect level.
The input file should have the structure:

    M1 = MASS ATOM 1
    M2 = MASS ATOM 2
    D0 = DISSOCIATION ENERGY IN HARTREE
    LEVELS = COMMA SEPARATED LIST AND HYPHENATED RANGE REQUESTING THESE VIBRATIONAL LEVELS. E.g. 1,2,3-5
    GUESS = COMMA SEPARATED LIST OF LEVEL:GUESS EIGENVALUES CORRESPONDING TO THE LEVELS ABOVE E.g. 1:-10,2:-8
    [
        SEPARATION (Bohr radii),TOTAL ENERGY (Hartree)
    ]

GUI

The GUI should be self explanatory, simply load the potential input file, modify the other data fields as required and hit “Solve the system” to run the solver. Text and graphical output can be saved using the respective buttons.

The GUI runs as expected for the first solution, but can become unstable when trying to load and solve subsequent potentials as the solver is not correctly reset between runs. As such it is recommended to restart the GUI for each system. If anyone finds the GUI routinely useful please let me know and I will fix this issue to allow sequential runs.

Solver as a python module

The Solver.py module can be imported into other scripts as with any other module. A Details object should be created by calling the readInput function providing path to the input file. The Details object returned from this function can be passed as an argument to the driver function which returns the solved system in a new Details object. Full documentation of the Details object and the properties containing the solution eigenvalues and wave functions can be found in the Solver.py file, and the example script.


License

This code is shared under the MIT liscence Copyright 2016 James Furness.
You are free to use, modify and distribute the code, though recognition of my effort is appreciated!