# Calculation of NMR shieldings using simple magnetic balance¶

In this tutorial we look at 4-component relativistic calculations of NMR shielding constants. For best results we recommend to use the scheme of simple magnetic balance (sMB) in conjunction with London orbitals [Olejniczak2012].

## Basic theory of magnetic balance¶

We encourage the reader to consult the original paper [Olejniczak2012] for full details; here we just provide some key points.

In the absence of magnetic fields the coupling between the large and small components of the Dirac equation is

The exact coupling can not be implemented *per se* at the basis set level due to
the energy dependence of the coupling. One therefore employs the
non-relativistic limit of this coupling

valid for positive-energy solutions (and extended nuclear charge
distributions). Construction of the small component basis in this manner
assures the correct representation of the kinetic energy operator in the
non-relativistic limit, and the coupling is therefore referred to as *kinetic
balance* [Stanton1984]. In a 2-spinor basis kinetic balance can be
implemented in such a manner that there is a 1:1 ratio between the sizes of the
large and small component basis; this is referred to as *restricted kinetic
balance* (RKB).

In a scalar basis the straightforward implementation of RKB is not possible and
one typically resorts to *unrestricted kinetic balance* (UKB)

basically generating the small component basis functions as derivative of the large component ones.

Magnetic fields are introduced through minimal substitution

which manifestly modifies the coupling between the large and the small components. The non-relativistic limit now reads

and is referred to as *magnetic balance* [Aucar1999]. Magnetic balance needs
to be taken into account when constructing the first-order correction to the
wave function upon inclusion of an external magnetic field and is not assured
in a RKB basis [Pecul2004]. In the cited paper [Olejniczak2012] we show that
magnetic balance can be obtained in a Gaussian basis combining UKB with London
orbitals. However, UKB generally leads to a small component basis larger than
the large component basis set from which it was generated and may increase
computational cost as well as introduce linear dependencies in the basis set as
compared to RKB. We therefore propose calculations of NMR shielding parameters
by a simple two-step procedure described in the next section.

## A sample calculation¶

We will consider a 4-component relativistic density functional calculation of
the NMR shielding constant of hydrogen fluoride using the KT2 functional
[Keal2003]. The molecular input file `hf.mol` reads:

```
DIRAC
C 2 A
1. 1
H 0.0000 0.0000 0.0000
LARGE BASIS aug-cc-pVTZ
9. 1
F 0.0000 0.0000 0.9168
LARGE BASIS aug-cc-pVTZ
FINISH
```

### Step 1: SCF calculation¶

The first step is the generation of the Kohn–Sham wave function in a RKB basis. Although DIRAC employs scalar Gaussian basis functions restricted kinetic balance is by default imposed when transforming to orthonormal basis [Visscher2000].

The menu file `scf.inp` reads:

```
**DIRAC
.WAVE FUNCTION
**INTEGRALS
*READIN
.UNCONTRACT
**HAMILTONIAN
.DFT
KT2
**WAVE FUNCTION
.SCF
*SCF
.CLOSED SHELL
10
**END OF
```

and we make sure to recover the coefficient file `DFCOEF` when running the calculation:

```
pam --inp=scf.inp --mol=hf.mol --outcmo
```

### Step 2: Response calculation¶

In the second step we start from the RKB coefficient obtained in the SCF calculation, but DIRAC will add the UKB complement
and
calculate the NMR shielding constant using unrestricted kinetic balance
combined with London orbitals. The menu file `response.inp` accordingly reads:

```
**DIRAC
.PROPERTIES
**GENERAL
.RKBIMP
**INTEGRALS
*READIN
.UNCONTRACT
**HAMILTONIAN
.URKBAL
.DFT
KT2
**WAVE FUNCTION
.SCF
*SCF
.CLOSED SHELL
10
**PROPERTIES
.SHIELDING
*NMR
.LONDON
**END OF
```

and the `pam` command is:

```
pam --inp=response.inp --mol=hf.mol --incmo
```

### Comparison with RKB calculations¶

With the above computational procedure we obtain the following NMR shielding constants:

```
isotropic shielding
----------------------------
atom total dia para
---------------------------------------------------
H 29.9502 21.8737 8.0765
F 415.8316 480.2015 -64.3699
---------------------------------------------------
```

By not adding the UKB complement (deleting .RKBIMP and .URKBAL), that is, using RKB, gives:

```
isotropic shielding
----------------------------
atom total dia para
--------------------------------------------------------
H 28.8028 20.7264 8.0765
F 323.2672 387.6459 -64.3787
--------------------------------------------------------
```

If we in addition switch off the use of London-orbitals (deleting .LONDON) we obtain:

```
isotropic shielding
----------------------------
atom total dia para
------------------------------------------------------
H 28.2920 20.2563 8.0357
F 322.9065 387.1375 -64.2310
------------------------------------------------------
```

The experimental value for the absolute isotropic shielding available for H in HF is 28.72 ppm [Chesnut1994].