Daniel J. Kowalski

Step 1 | Setup for Self-Consistent Field (SCF) Calculation

This section outlines the structure of an input file for Quantum ESPRESSO's plane wave executable (pw.x). Once written, several parameter values can then be optimised to ensure an appropriate level of accuracy in the simulation (see Section 2).

Documentation for the pw.x executable is available from the Quantum ESPRESSO website.

Cell parameters were sourced from experimental crystal structures where possible with the Materials Project as a fallback. For rutile, the experimental structure was sourced from the Inorganic Crystal Structure Database (ICSD 142924, Deposition #: 2187406), originally recorded in Poschmann et al. (DOI: 10.1002/ejic.202000555) The CIF file was interpreted using the Gemmi, pymatgen, and Atomic Simulation Environment (ASE) packages for Python, with the cell dimensions and atomic positions used directly for the SCF calculation.

Pseudopotentials were taken from the GBRV high-throughput pseudopotentials set produced by the Vanderbilt lab. The relevant reference is Comput. Mater. Sci., 2014, 81, pp.446 (PPs available from: https://www.physics.rutgers.edu/gbrv/).

The K point grid was chosen to provide a reasonably dense sampling, as this was possible for rutile's small unit cell without being too computationally expensive. As the unit cell is tetragonal, and the c-axis is shorter than the a- and b-axes, a finer sampling was used to maintain comparable sampling density in the reciprocal lattice.

Step 2 | Optimising Key Parameters: ecutwfc, ecutrho, degauss, K points

In order to optimise the calculation, a series of Linux shell scripts were used to vary the parameter of interest with a 'for do' loop. Once all calculations for a given parameter were run, the total energies were compared. Convergence was assumed once the variation from one value to the next fell below 0.005 Ry (with any remaining variation assumed to be due to temperature). Initially, the 'ecutwfc' value was optimised, then 'ecutrho', and finally the K points grid. For metallic systems, the 'degauss' parameter was also optimised simultanteously with the K points grid.

For rutile, 'ecutwfc' converged at 50 Ry. 'ecutrho' converged at 200 Ry. This was slightly unexpected as the Ti PP used is an ultrasoft potential, which often require higher ecutrho values than the typical 4 x ecutwfc. K point grid converged at 4 4 6 0 0 0 (corresponding to ./opt/rutile_kpts4.out in the example above).

Step 3 | Variable Cell Relaxation

Finally, a vc-relax calculation was run using the optimised values from step 2 to optimise the cell dimensions and atomic positions. This calculation required several additional parameters and the includion of '&IONS' and '&CELL' sections. In the main, the additional parameters provide convergence thresholds for the calculation. Many of these were added explicitly when earlier attempts failed to converge for rutile, but may not be necessary for other calculations. These are set to the pw.x documentation's suggested values, and while they are made explicit here, may have been included implicitly in earlier calculations.

The '&IONS' section here explicitly sets the algorithm used during relaxation ('ion_dynamics') and drops the intial ionic displacement used ('trust_radius_ini'). This latter change was made when earlier attempts failed to converge.

The '&CELL' section here explicitly sets the algorithm used during relaxation ('cell_dynamics'), pressure for the calculation ('press') and pressure convergence threshold ('press_conv_thr'). The degrees of freedom parameter ('cell_dofree = 'ibrav'') requires the calculation to retain the tetragonal lattice. 'cell_factor = 2.0' is a required parameter for construction of pseudopotential tables while optimising cell dimensions.

Results | Rutile

To assess the efficacy of my relaxation of rutile, the cell parameters were compared against experimental data and materials project simulations. In general, the unit cell dimensions from my calculation are slightly closer to the experimental data than those proposed by the materials project. The materials project calculated values for atomic positions are in slightly better agreement with the crystallographic data than mine are, however the difference is slight. My values are within about 1% of both sets of values for all measurements.

Precision used for each parameter is set equivalent to the least precise measurement. Atomic positions are given in fractional coordinates.