How to make a k-point convergence test for a bulk structure with a lower k-point in the z-direction?

Question

Let's say I want to make k-point convergence test for graphite. And let's say it converges at 12x12x4.

Do I first need to do a k-point convergence test for equal k-points on xyz directions then make another k-point convergence test with changing k-point only along z-direction?

I am asking because I don't want to use 12x12x12 k-points if it already converges for 12x12x4.

Answer 1

1. If you have already obtained satisfactory convergence with a (relatively) sparse k-point grid, there is no motivation to go for a denser grid. So if you have already achieved convergence with 12x12x4, there's no need to go to 12x12x12.

2. If you are talking about graphene, which is 2-D, there is no need to sample points along the out-of-plane direction.

3. There are some more nuances, which you haven't mentioned in the question, but I will comment anyway. Even grids usually have better convergence than odd grids because you avoid sampling the high-symmetry points. So you should typically try to use even k-meshes for your BZ integrations. But there is an exception (for materials with hexagonal BZs) in the case of graphene, you should not use a gamma-centered (even) grid because you tend to generate points outside the brillouin zone.

Answer 2

I agree with the answer provided by Xivi76. I just wanted to add that some codes have a very nice functionality that facilitates convergence: rather than explicitly writing out a $\mathbf{k}$-point grid $n_1\times n_2\times n_3$, in which in principle you have to converge three values; you can instead specify a $\mathbf{k}$-point spacing or density, in which case you only have to vary one value. This is useful in the following situations:

1. The cell lattice parameters are of different length, like in your graphite example. In this case, specifying the spacing directly generates fewer $\mathbf{k}$-points along Brillouin zone directions that are shorter (corresponding to longer real space lattice parameters).
2. You are interested in performing a series of supercell calculations (e.g. to study phonons). In this case, specifying the spacing directly allows you to construct consistent $\mathbf{k}$-point grids for the various supercell sizes and shapes. This consistency is not perfect of course, because it may not be possible to exactly divide the dimensions of a $\mathbf{k}$-point grid precisely into a smaller Brillouin zone arising from a given supercell, but it does the best possible job.

As examples, Castep uses the keyword KPOINT_MP_SPACING, while Vasp uses the keyword KSPACING.

Answer 3

• As far as my knowledge this a relationship between k-points and lattice constants values. I'll give you an example of a layered hexagonal material WS2:

Lattice constants : a=3.17 b=3.17 c=12.41 ; so c/a = 12.41/3.17=3.91

K-points : If I chose kx=12, ky will be equal to 12 but kz should be equal to an integer close to the kx divided by c/a. I mean kz=12/3.91=3.06. So the consistent k-points will be 12x12x3

Let's say that I chose kx=9, ky=9 and kz will be =9/(c/a)=9/4=2.25. So we'll have 9x9x2

• Remember it is linked with your c/a ration and k-points should be integers.
• You have to conduct first K-points convergence tests for multiple values (e.g.: 9x9x2, 10x10x3,11x11x3 ,12x12x3 ...etc) and check where the energy is minimum.

Note : for cubic systems the k-points are equal, because lattice constants are equal.

Answer 4

I will try to give the most practical answer, the reality of "is this converged" is that you cannot know without checking by going past it. You say that you would like to save time doing these calculations, but the worst loss of time is sometimes the loss of your own human time as you get confusing results later.

I am unsure how you came to the number of 4 (cell vector * kpts?), but lets see using GPAW if that makes sense. I am using what are likely unconverged settings with an unrelaxed structure. Even then, lets see if we see a convergence on a (12, 12, X) kpoint set.

As you can see, there is something strange going on around 6 or 7 kpts in the z direction but otherwise looks converged at 4. Tightening my settings might remove that strange bump. If I saw this and could not correct it, I would probably choose to run initial optimizations at (12, 12, 4) and final optimizations at (12, 12, 8).

I strongly recommend you do the same and check what you see. If you are unsure what convergence looks like you could provide us with a similar graph. If you have a property of the system you are investigating, use that as well as total energy.

GPAW version 20.1.0 and ASE version 3.19.2 were used to generate this data. Here is the script I used.

from ase.build import bulk
from ase.visualize import view
from gpaw import GPAW, PW
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
import numpy as np

a = 1.42
c = 6.70

kpt_min = 1
kpt_max = 16

atoms = bulk("C", crystalstructure="hcp", a=a, c=c)

x = np.arange(kpt_min, kpt_max)
e = np.zeros((kpt_max-kpt_min))

for index, kpt in enumerate(x):
    calc = GPAW(mode=PW(350), kpts=(12, 12, kpt), occupations={'name': 'fermi-dirac', 'width': 0.05})
    atoms.calc = calc

    e[index] = atoms.get_potential_energy()

plt.plot(x, e, linestyle="-", marker="o")
plt.show()