Calculating and using Dynamical Matrix

|U2022.12|

Downloads & Links

DynamicalMatrix
Basic QuantumATK Tutorial
ATK Reference Manual

Introduction

The dynamical matrix is an object that is used to calculate many properties related to dynamics of the crystal lattice e.g Mobility, RamanSpectrum, DeformationPotential, InelasticIVCharacteristics and more. In this tutorial you will learn how to calculate the dynamical matrix and use it to perform a simple phonon bandstructure calculation. In particular, you will learn how to set up the calculation and which parameters are critical. You will also learn how to optimize the performance of this type of calculations by an optimal parallelization strategy.

• You will first perform am DFT-LCAO calculation with the LDA exchange-correlation potential. You will go through all the relevant parameters used in a phonon calculation.
• Towards the end, you will use the classical potentials available in the ForceFieldCalculator and compare the results.

Create the Workflow

Go to the Builder and import the Silicon (alpha) structure from the built-in database. Then send the structure to the Workflow Builder and add the following script blocks:

1. LCAOCalculator
2. OptimizeGeometry
3. DynamicalMatrix
4. PhononBandstructure

In the end the workflow should look like this

Now we will go through the settings of each block and edit the important parameters.

LCAOCalculator Settings

To get good results for the dynamics one often has to converge with respect to a few key parameters:

• Open the LCAOCalculator block.
• In the Main section, in Basic Parameters select LDA under Family, FHI under Pseudopotential and SingleZetaPolarized under Basis Set. This is a minimal basis set that is used to speed up the calculations.
• In the Main section, in Numerical Accuracy set Broadening to 300 Kelvin, Density mesh cutoff to 50 Hartree and k-point density to 6 Angstrom.

• In the Iteration Control section, set a Damping factor of 0.5. Also reduce the Tolerance to 0.000001.

Tip

• It can often be necessary to increase mesh cutoff and k-point density for DynamicalMatrix calculations to ensure convergence of the phonon bandstructure.
• For non-metals in general we recommend using 300 Kelvin broadening.
• For DynamicalMatrix calculations it is generally recommended to use a tighter tolerance than the default one.
• It is possible to reduce the number of iteration steps by tuning the damping factor and number of history steps used.

Lattice optimization Settings

Open the OptimizeGeometry block, and modify the parameters as in the figure below to perform a force and stress minimization.

Tip

• We generally recommend reducing the target stress and forces when calculating the dynamical matrix.
• To get good structures for calculations of electronic properties the default settings are usually sufficient.

Dynamical matrix Settings

The dynamical matrix is the second derivative of the energy, corresponding to the first derivative of the forces. The first derivative of the forces are calculated using a finite difference scheme, where each symmetrically unique atom in the central cell is displaced along all cartesian directions, also called frozen phonon calculations.

For crystals with small unit cells, a periodically repeated super-cell of the unit cell is often needed to accurately account for long-range interactions in the dynamical matrix. QuantumATK can handle this automatically. By default QuantumATK uses enough repetitions so that each atom is at least 4 covalent bond radii away from its translated copies. This normally includes atoms that are up to 4 bonds away. This is more than enough for most systems. It is also possible to directly set the number of repetitions used in each direction. This can especially be important in DFT calculations, which can become quite heavy depending on the number of repetitions. You can see more details in notes of the manual entry for DynamicalMatrix.

• In the Dynamical Matrix object set Repetitions to Custom. Set n_A, n_B and n_C to 7. Press OK

• Open the PhononBandstructure block and choose 400 points per segment and the G, X, K, G, L Brillouin zone route.

To test the effect of using different numbers of repetitions, we repeat the calculation setting Repetitions to 3x3x3 and 5x5x5. The full workflow can be downloaded here Silicon_phonon_bandstructure.hdf5.

Running the calculation

You are now ready to run your DFT phonon bandstructure calculation. This can be done by sending the script to the Job Manager using the Send To icon .

However, there are a couple of details to consider before running the job:

Parallelization

In QuantumATK, the dynamical matrix can be calculated by parallelizing over displacements. The number of displacements needed depends on the system and can be calculated by opening the Dynamical object and pressing the Calculate… button under Number of Displacements. This means that the optimal performance can be obtained by using a number of MPI processes equal to the number of displacements + 1 process reserved for delegation. This can be controlled using the Processes per displacement.

For this system, six different displacements need to be calculated.

System size

The three systems run in this tutorial, 3x3x3, 5x5x5 and 7x7x7 repetitions, consists of 54, 250 and 686 Silicon atoms, respectively. Besides the running time (see below) you should also consider the amount of memory required. A bulk configuration with 686 silicon atoms run with LDA and with a SZP basis set will take about 3 GB of memory. Each displacement, i.e. each MPI process, uses the same amount of memory. This means that in this case (6 MPI processes) if you run on a single machine you will need at least 18 GB of memory!

Timing

The table below reports the total time required to run a phonon bandstructure calculation with the parameters set above. The calculations are run in parallel with 20 MPI processes:

Repetition	3x3x3	5x5x5	7x7x7
Total time	1m29.82s	22m26.95s	1h39m44.40s

Analyzing the results

When the job is done you will find in the Data view the PhononBandStructure object. Highlight it and press Open to display it.

The figure below reports results for the 3x3x3, 5x5x5 and 7x7x7 systems defined above.

You immediately notice that the 3x3x3 repetition is too small to give even qualitatively good results, while the intermediate size system is already a pretty decent approximation. However, you can also notice that a small difference exists between the acoustic modes near the X-point for the 5x5x5 and the 7x7x7 systems.

The dynamical matrix is now stored together with the BulkConfiguration. This means you can read the dynamical matrix and perform another PhononBandstructure calculation (on a different Brillouin path for example) or a PhononDensityOfStates calculation in just a few seconds.

Speeding up the calculation with ForceFields

QuantumATK offers a large set of classical potentials in the ForceFieldCalculator. These potentials are well suited to study vibrational properties of materials, such as the phonon bandstructure of silicon calculated in this tutorial.

Note

Classical potentials are available only for a limited amount of systems. As your system becomes more complex including e.g. interfaces, defects or many elements you will often find that no relevant classical potentials exist. In this case QuantumATK offers the ability to train your own machine learned potential using Moment Tensor Potentials. This comes as part of the Elite license package along with many user friendly tools for automated generation of training structures.

Send the silicon bulk configuration to the Workflow Builder. Add a ForceFieldCalculator. Now a stress optimization and phonon bandstructure calculation can be set up as explained above.

The calculation will take just a few seconds on a normal desktop machine.

In ForceFieldCalculator select the most appropriate potential available for silicon, such as the Tersoff or the Stillinger and Weber potentials.

Each potential is provided with a literature reference, which you should check to determine if the potential is well suited for your system and for the property you are interested in.

Set up the OptimizeGeometry block and the PhononBandstructure analysis object exactly as for the LCAOCalculator case and run the calculation.

The phonon bandstructure of silicon calculated with the Tersoff_si_2005 potential [2] is reported in the figure below.

References

[1]	Giannozzi et al., “Ab initio calculation of phonon dispersions in semiconductors”, Phys. Rev. B 43, 7231 (1991).

[2]	F. H. Stillinger and T. A. Weber, “Computer simulation of local order in condensed phases of silicon”, Phys. Rev. B, 31, 5262 (1985).