In this example, we will set up and run a driver calculation in Orca. In this type of calculation, one variable (an angle or bond length) is varied while the energy is minimised. This type of calculation can be useful for exploring reaction pathways. In this example, we will use the dissociation of a carbonyl ligand from a rhodium-phosphine complex. These types of complexes are common in industrial chemistry for a variety of catalytic reactions; however, typically more complicated phosphine or phosphite ligands, either mono- or bidentate, are used. For computational simplicity, I am using phosphine (PH3) as the ligand.

Setting up the calculation

For simplicity, we will use a split valence basis on all atoms except rhodium. On rhodium, we will use a Stuttgart-Dresden ECP with the def2-tzvp basis set modelling the valence electrons. The calculation uses BP86 and Grimme’s D3 dispersion correction with Becke-Johnson damping. Below is the input for this driver calculation.

! noautostart slowconv tightscf d3bj opt bp86 pal2 def2-svp def2-svp/j ri ecp{def2-sd=[Rh], def2-tzvp, def2-tzvp/j}
%geom scan
B 0 11 = 1.909, 4.0, 8
*xyz 0 1

Concerning the topic of this post, i.e., the driver calculation, the most important part of this input file is the part describing the bond length that we are ‘driving’. In this example, the Rh-C bond length. In the XYZ coordinates, Rh and C (of carbonyl) are the first and twelfth atoms, respectively. Orca numbers atoms from 0 (zero) upwards, so it is important to remember that the nth numbered atom is the n-1th atom in Orca. Here, the command ‘B 0 11 = 1.909, 4.0, 8’ indicates that a bond (B) is changed from 1.909 to 4 Å in eight steps, i.e., 2.6 Å per step.

In addition, the ‘simple’ command line contains several additional elements. The first command ‘noautostart’ indicates that we do not want to read from an existing wavefunction file, should it exist. ‘slowconv’ and ‘tightscf’ are commands that control the SCF calculation, and these are useful for transition metal complexes. The command ‘d3bj’ indicates that we want to model dispersion interactions using Grimme’s D3 correction with Becke-Johnson damping, and ‘opt’ indicates that this is an optimisation, that is, optimising everything other than the coordinate being driven. ‘pal2’ indicates that this is a parallel calculation, running on two openmpi processes. ‘def2-svp def2-svp/j ri ecp{def2-sd=[Rh], def2-tzvp, def2-tzvp/j}’ are commands related to the basis sets on the atoms in this system. A split basis (def2-svp) is being used on all atoms; however, on Rhodium, an electron core potential is being used to model the core electrons, and a triple-zeta basis is used for the valence electrons. In addition, the RI approximation is used to speed up the calculations, and this si the reason for the basis sets appended with ‘/j’, i.e., these are the density fitting basis sets.

About the calculation

The GIF below shows the eight steps of the relaxed driver calculation; for fun, it is a stereogram. The Rh-C bond was lengthened in steps of about 0.3 A. At each step, an optimization was carried out, but the Rh-C distance was kept constant.


As the carbonyl ligand moves away from the rhodium atom, the Rh-C-O bond angle changes. Starting off linear at 177.2 degrees; it becomes bent with an angle of about 109.

Charting the Rh-C-O bond angle versus the Rh-C bond length.

Of course, after a certain point there is no interaction between Rh and the C of the carbonyl ligand; this is given by the dissociation limit of the complex. The Rh-C-O tilt angle shows how the bonding interaction changes from a ligand sigma bonded with π-backbonding to a free ligand, and this is nicely shown by the shortening of the C-O bond length from 1.16 A to 1.14 A as the carbonyl ligand dissociates.

Analysing the data

Energy profile for the driver calculation.

The results of the calculation are shown to the left. As you can see, the CO-metal bond elongates, followed by a tilting of CO and its slow drift away from the complex. The energy of the complex increases until it reaches the dissociation limit; here, the energy levels off to a constant value.

In fact, we can fit the energy profile to the Morse potential and get an idea of the dissociation limit and the equilibrium distance for the Rh–C bond length.

The Morse equation is y = D \times(1-exp(-a\times (x-e)))^2 where D is the dissociation limit, a is the well depth, and e is the equilibrium distance.

Fitting of a Morse curve to the driver data.
Fitting of a Morse curve to the driver data.

Fitting this equation to the data, The morse curve gives values for D, a, and e of 26.92 kcal/mol, 2.81 kcal/mol, and 1.92 A, respectively.


Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s