I want to estimate the binding energies of various molecules. What are the analytical and numerical tools I can use for this purpose? In particular I'm interested in bonding energies of dyes to substrates (example).

  • 3
    $\begingroup$ Are you only referring to pen and paper methods, or are you also interested in experimental procedures? It would probably be helpful if you could give some examples of the molecules (and bonds) you are looking into to estimate the size and what methods are suitable. $\endgroup$ Commented Jan 29, 2015 at 2:23
  • 1
    $\begingroup$ This will depend heavily on the accuracy you need and the size of the molecules. Binding one protein to another is very different than dealing with the benzene dimer. $\endgroup$ Commented Jan 29, 2015 at 2:35
  • $\begingroup$ @Martin, I've updated the question. $\endgroup$
    – Sparkler
    Commented Jan 29, 2015 at 2:45
  • 1
    $\begingroup$ @Jori I consider computational chemistry to be a pen and paper method, at least in principle. The computer just does the job. While Gaussian might be suitable for a lot of things, i believe for the OP's intention it might not be the best fit. I might try to add an answer tomorrow, but I welcome anyone with more insight to be faster. $\endgroup$ Commented Jan 29, 2015 at 16:33
  • 1
    $\begingroup$ @Martin How were the results of Gaussian (the very first version) ever verified? Prof. John Pople simply got himself an army of grad students who worked overtime crunching numbers with paper, pencil and those old cumbersome calculators that had those huge levers you had to 'sha-shunk!' every time you wanted to perform an operation. I'm so glad I wasn't practicing comp. chem. during those times. ;) $\endgroup$ Commented Jan 30, 2015 at 15:45

1 Answer 1


The steps for predicting a binding energy is relatively straightforward in computational chemistry but the process itself can become extremely difficult depending on

  • The types of molecules involved (size, complexity, etc.)
  • The desired accuracy of your final answer
  • The environment that the molecules reside in (gas phase, aqueous phase, etc.)

For very small molecules in the gas phase, the process is easy and routine assuming that electron correlation is trivial. Aqueous phase environments, complex/large molecules and high accuracy quickly makes the problem more difficult if not impossible.

I have written up a previous answer which addresses how to predict an interaction energy ($E_{\mathrm{int}}$) which is completely different from a binding energy ($E_{\mathrm{bind}}$). However, much of the concepts will remain the same. I will try to outline here what in fact is the difference between an interaction energy and a binding energy. Unfortunately, 'binding energy' is tossed around a lot in the literature and consequently grossly misused. I hope to clarify the subtle differences between these two terms in my answer to your question.

I will focus my discussion around an elementary system of two water molecules to make this conceptually easier. In the gas phase, an isolated water molecule will adopt a particular geometry (at 0K and in a vacuum). We can determine what this minimum energy configuration will be at a particular level of theory by approximating the Schrodinger equation. The geometry that gives rise to the lowest possible energy for an isolated water molecule will be the geometry that you will observe under those conditions. Remember that the internal coordinates of water include things like the O-H bonds, and the H-O-H angle. We don't need to guess at what these parameters will be in the end, solving the Schrodinger equation will determine that.

Now, we know exactly what an isolated water molecule looks like and what the corresponding energy of that water molecule will be. Both of these components are equally important when it comes to determining a binding energy because (as the name suggests) we are trying to determine an energy. So the energies that we start off with are crucially important.

We can now consider what happens when two water molecules interact (form a hydrogen bond). What will be the binding energy of these two water molecules (i.e. what will be the strength of this bond that forms between then?). It is common to approximate the strength of this bond by simply computing an interaction energy (its computationally easier to determine $E_{\mathrm{int}}$ rather than $E_{\mathrm{bind}}$). How would we compute an interaction energy? In this case, you wouldn't need to figure out what geometry or energy an isolated water molecule in its lowest energy configuration. You would simply need to determine what the geometry of two interacting water molecules (i.e. dimer) would be (and corresponding energy) as well as the energy of each water molecule isolated from each other in the geometry that it adopts in the dimer.

I'll explain. We can determine the lowest energy configuration of a single water molecule in the gas phase. This geometry is inevitably different than a water that is not isolated... a water that is interacting with something else such as another water. If you want to get an interaction energy, you only need...

  • Energy of the lowest energy configuration of the water dimer
  • Energy of each water (in the water dimer) isolated from each other but only at the geometry that the water takes on within the dimer.

In this case, you are NOT allowing each isolated water molecule to 'relax' into its lowest energy configuration. You are simply going to neglect the contributions of this energy relaxation effect to the overall bond strength. Therefore,

$E_{\text{int}} = E_{\text{opt-dimer}} - (E_{\text{mon1}} + E_{\text{mon2}})$

where mon1 is the first water molecule in your dimer and mon2 is the second water molecule in your dimer. (opt-dimer means optimized dimer)

However, in the case of a binding energy, you MUST account for the relaxation effects that occur when you break the dimer apart and allow the waters to become isolated from each other (and everything else). To do this, you simply carry out what I've stated above... you simply optimize each water and then take its corresponding energy. (Determining an interaction energy means you don't have to do this extra optimization(s)... you would simply do a single point energy on each monomer). Clearly for a water dimer, you would only do one optimization as each water will become identical when isolated. For a heterodimer, you would have to optimize each piece separately.

To determine the binding energy, you would use the following equation:

$E_{\text{bind}} = E_{\text{opt-dimer}} - (E_{\text{opt-mon1}} + E_{\text{opt-mon2}})$

Of course the problem becomes more complex when taking into account basis-set superposition error (you will need to perform counterpoise corrections though there exists a heated debate in the literature about the usefulness of this type of correction). Also, zero-point vibrational energy corrected binding energies will require that the zero-point energy of the dimer and each monomer will need to be determined (which comes from taking the 2nd derivative of the wavefunction which is computationally demanding).

As far as tools go, you can choose from a number of available quantum software packages. Gaussian 09 is a popular choice but very expensive. Q-Chem (paid) has the benefit of having the most beautiful, gorgeous, well-laid out documentation of anything I've ever seen (makes for a good supplemental textbook to boot). NWChem is opensource. Google can help you find plenty of options.

A quick note: Basis set superposition error (BSSE) should NOT be referred to as an error as 'error' suggests deviation from the RIGHT answer. BSSE has been shown to actually make non-corrected answers even WORSE.

One last important note: When determining differences of energies (such as you do to get an interaction or binding energy), ALL energies being compared must be evaluated at the SAME LEVEL OF THEORY (method/basis set, etc.).

  • $\begingroup$ I'm looking at NWChem's sample input files. If I want to compare "binding energies" of various combinations, can I use the SCF method on two molecules separately and substract from the SCF value for the dimer? $\endgroup$
    – Sparkler
    Commented Jan 30, 2015 at 16:25
  • 1
    $\begingroup$ Yes, as long as you are comparing energies that were all determined at the same level of theory (method/basis etc.). Be sure to convert Hartrees to kJ/mol (or kcal/mol if you're still stuck in the 1970s). Hartrees is unnecessary and inefficient when taking differences of energies (Hartrees are appropriate for absolute energies but absolute energies are completely inaccurate anyway). $\endgroup$ Commented Jan 30, 2015 at 16:29
  • $\begingroup$ For what it's worth, you may consider using SAPT to compute interaction/binding energies directly. It also allows for a physical decomposition of the nature of the interaction. PSI4 contains it, is free, and is apparently pretty fast. sirius.chem.vt.edu/psi4manual/4.0b5/sapt.html $\endgroup$
    – jjgoings
    Commented Jan 30, 2015 at 16:38
  • $\begingroup$ @LordStryker ...I like kcal/mol :/ $\endgroup$
    – jjgoings
    Commented Jan 30, 2015 at 16:38
  • $\begingroup$ @jjgoings Cautionary note. The SAPT routines for anything greater than a dimer interaction is broken and shouldn't be used. 3-body and above methods are utterly wrong (last I checked). $\endgroup$ Commented Jan 30, 2015 at 16:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.