Silicon phonon bandstructure

Introduction

In this tutorial you will learn how 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 ATK-DFT 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 ATK-ForceField engine and compare the results.

Lattice optimization

The first, and very important, step is to perform a stress optimization of your system. This removes negative frequencies in the phonon spectrum, which are in general due to a non-fully relaxed structure.

  1. Go to the Builder builder_icon and import the Silicon (alpha) structure from the Database database_icon.
  2. Send the structure to the Script Generator script_generator_icon and add a New Calculator:
  • Select ATK-DFT and LDA (default options).
  • Use a 13x13x13 k-point sampling.
  • In Iteration control parameters, use a damping factor of 0.5. This value reduces the number of iterations needed to converge the SCF state for this particular system.
  • Use a SingleZetaPolarized basis set to speed up the calculations.
  • Add an Optimization ‣ OptimizeGeometry block, and modify the parameters as in the figure below to perform a force and stress minimization.
  • Change the output filename to Silicon_opt.nc, send the script to the Job Manager job_manager_icon, and run it.
../../_images/Silicon-optimizegeometry_stress.png

The calculation will take just a few seconds. Afterwards, you will find the optimized BulkConfiguration in the LabFloor with the object id glD001. Select this object and click on Text representation to get the optimized lattice parameters, 5.4787 Angstrom.

../../_images/Silicon-LabFloor_optimized_structure.png

Set up the phonon bandstructure calculation

Drag and drop the optimized BulkConfiguration to the Scritper script_generator_icon.

  1. Double click on New Calculator and modify the k-point sampling to 1x1x1.
    • This is the k-point sampling that will be used for the calculation of the dynamical matrix, which is performed with a repetition of your configuration (see below). See also K-points used for ATK-DFT and ATK-SE.
  2. Add a PhononBandStructure analysis object.
    • Choose 400 points per segment and the G, X, K, G, L Brillouin zone route.
    • Change the output filename to SiliconPhBS.nc.
    • In order to access to the dynamical matrix parameters (next step), set Script detail to “Show defaults” and send the script to the Editor editor_icon:
../../_images/SiliconPhononBS_send_to_editor1.png
  1. Phonons are calculated from the dynamical matrix of the system. 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 the system is displaced along each degree of freedom in the system, also called frozen phonon calculations.
    • The details of the calculation of the dynamical matrix is controlled by the DynamicalMatrixParameters object, which can be given as an option to the calculator. In order to take into account the dynamical matrix entries with neighbour cells, the system is repeated, controlled by the repeat parameter. The default value is “automatic”; ATK will try to determine the relevant repetition and this is reported in the output file.
    • You have access to the DynamicalMatrixParameters if you inspect the script with the Editor editor_icon, showing all the default parameters, as you did in the previous step.
    • Locate the dynamical_matrix_parameters section in your script. The default parameters are:
dynamical_matrix_parameters = DynamicalMatrixParameters(
    atomic_displacement=0.01*Angstrom,
    repeats=Automatic,
    acoustic_sum_rule=False,
    symmetrize=True,
    )
  1. Perform a few modifications:
    • acoustic_sum_rule: Change it to True in order to invoke the acoustic sum rule.
    • repeats: With this parameter you select how many times your system is repeated in the x, y and z direction for the calculation of the dynamical matrix.
      • The default value is Automatic and this is suitable for most systems. However, you should pay attention to the resulted repeated structure, since the calculation can become pretty heavy if want to use DFT.
      • In this tutorial, you will calculate the phonon bandstructure for a 3x3x3, 5x5x5 and 7x7x7 repetitions with DFT.

For the 5x5x5 repetition, your dynamical_matrix_parameters will now look like this:

dynamical_matrix_parameters = DynamicalMatrixParameters(
    atomic_displacement=0.01*Angstrom,
    repeats=(5,5,5),
    acoustic_sum_rule=True,
    symmetrize=True,
    )

Running the phonon bandstructure calculation

Save the script. You are now ready to run your DFT phonon bandstructure calculation.

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

Parallelization

In ATK, the dynamical matrix is calculated by parallelizing over displacements. The number of displacements is equal to three times the number of atoms (one displacement per each x, y, z directions). For each displacement, two calculations are performed sequentially, corresponding to + and – directions. This means that the optimal performance is obtained by using the same number of MPI processes as the number of displacements.

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

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 6 MPI processes on a Xeon e5472 3.0GHz machine:

Repetition 3x3x3 5x5x5 7x7x7
Total time 6m21.44s 41m16.41s 3h06m53.15s

Analyzing the results

When the job is done you will find in the LabFloor the PhononBandStructure object. Select it and use the Bandstructure Analyzer tool to display it.

../../_images/Silicon_LabFloor_phononbandstrcture.png

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

../../_images/Silicon_phononbandstructure_repetitrions.png

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 fine approximation. However, you can also notice that a small difference exists between the acoustic modes near the \(\Gamma\)-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:

../../_images/Silicon_PhononBandstructure_Silicon-analysis_from_file.png

The following figures show the phonon bandstructure and DOS, which you can compare to the results reported in Ref. [1].

../../_images/Silicon_PhononBandstructure_Silicon-compare_phononbandstructure.png ../../_images/Silicon_PhononBandstructure_Silicon-compare_densityofstates.png

ATK-ForceField

In ATK 2014 and newer, a large set of classical potentials is included in the ATK-ForceField engine. These potentials are well suited to study vibrational properties of materials, such as the phonon bandstructure of silicon calculated in this tutorial.

Send the silicon bulk configuration to the Scripter script_generator_icon, and set up a stress optimization and phonon bandstructure calculation as explained above, both in the same script.

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

../../_images/Silicon_PhononBandstructure_Silicon-atkclassical_scripter.png

In New Calculator, choose the ATK-ForceField engine and use the most appropriate potential available for silicon, such as the Tersoff or the Stillinger and Weber potentials.

Note

For ATK-versions older than 2017, the ATK-ForceField calculator can be found under the name ATK-Classical.

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.

../../_images/Silicon_PhononBandstructure_Silicon-atkclassical_scripter_SW.png

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

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

../../_images/Silicon_PhononBandstructure_Silicon-atkclassical_phononbandstructure.png

References

[1]
  1. Giannozzi et al., “Ab initio calculation of phonon dispersions in semiconductors”, Phys. Rev. B 43, 7231 (1991).
[2]
    1. Stillinger and T. A. Weber, “Computer simulation of local order in condensed phases of silicon”, Phys. Rev. B, 31, 5262 (1985).