Force constant potential¶

The force constant potential (FCP) is obtained through a systematic expansion of the energy in displacements as described [Brorsson2021].

Potential form¶

In the force constant potential, the potential energy is calculated as a Taylor expansion in terms of the atomic displacements \(\{u_i^{a}\}\) relative to a set of reference (equilibrium) positions as

\[U = U_2 + U_3 + U_4 + U_5 + U_6 + \cdots\]

with

\[\begin{split}U_2 &= \frac{1}{2!} \sum_{ij} \sum_{ab} \Phi_{ij}^{ab} u_{i}^{a} u_{j}^{b} \\ U_3 &= \frac{1}{3!} \sum_{ijk} \sum_{abc} \Phi_{ijk}^{abc} u_{i}^{a} u_{j}^{b} u_{k}^{c} \\ U_4 &= \frac{1}{4!} \sum_{ijkl} \sum_{abcd} \Phi_{ijkl}^{abcd} u_{i}^{a} u_{j}^{b} u_{k}^{c} u_{l}^{d} \\ U_5 &= \frac{1}{5!} \sum_{ijklm} \sum_{abcde} \Phi_{ijklm}^{abcde} u_{i}^{a} u_{j}^{b} u_{k}^{c} u_{l}^{d} u_{m}^{e} \\ U_6 &= \frac{1}{6!} \sum_{ijklmn} \sum_{abcdef} \Phi_{ijklmn}^{abcdef} u_{i}^{a} u_{j}^{b} u_{k}^{c} u_{l}^{d} u_{m}^{e} u_{n}^{f}.\end{split}\]

Here, \(\Phi_{ij}^{ab}\), \(\Phi_{ijk}^{abc}\), \(\Phi_{ijkl}^{abcd}\), \(\Phi_{ijklm}^{abcde}\), and \(\Phi_{ijklmn}^{abcdef}\) are the second-order, third-order, fourth-order, fifth-order, and sixth-order force constants. The indices \(i\), \(j\), \(k\), \(l\), \(m\), and \(n\) refer to the atoms and can take integer values from 0 to \(N-1\), where \(N\) is the number of atoms in the system. The indices \(a\), \(b\), \(c\), \(d\), \(e\), and \(f\) refer to the axes in the Cartesian coordinate system and can take integer values 0, 1, and 2, which correspond to the \(x\), \(y\), and \(z\) axes, respectively. In GPUMD, we only consider force constants up to sixth order.

File format¶

One needs to prepare several files related to this potential, which can be conveniently generated using the hiphive package [Eriksson2019], except for the driver potential file below (which is very easy to prepare).

driver file¶

The driver potential file for this potential model reads:

fcp number_of_atom_types <list of atom symbols>
highest_force_order highest_heat_current_order
path_to_force_constant_files

fcp is the name of this potential and tells the code that we are using the force constant potential.
number_of_atom_types is the number of atom types defined in the simulation model file.
<list of atom symbols> is a list of all the elements in the potential (can be in any order).
highest_force_order is the highest order of the force constants used in the potential. For example, when highest_order is 4, second-order, third-order, and fourth-order forces constants will be used (and should be prepared).
highest_heat_current_order is the highest order of the force constants used for heat current calculations. It can only be 2 or 3, which means using the second order force constants only or using both the second and third order force constants for heat current calculations, respectively.
path_to_force_constant_files is the path to the force constant files (see below). Important: There should be no trailing slash (/) after the folder name.

force constant files¶

The force constant data should be prepared in some files named:

clusters_order2.in
clusters_order3.in
clusters_order4.in
clusters_order5.in
clusters_order6.in
fcs_order2.in
fcs_order3.in
fcs_order4.in
fcs_order5.in
fcs_order6.in

These files should be in the folder you specified in the driver potential file (see above). If you only consider force constants up to the 4th order, you do not need the files with numbers 5 and 6. These files can be generated by the hiphive package. Here, we therefore do not describe the formats of these files.

equilibrium position file¶

Because this potential is defined in terms of the atomic displacements, one has to define the equilibrium (reference) positions of the atoms in the system. A file called r0.in is used for this purpose. This file should be in the folder you specified in the driver potential file (see above). The format of this file is:

x_0 y_0 z_0
x_1 y_1 z_1
x_2 y_2 z_2
x_3 y_3 z_3
...

where each line gives the position of one atom. The order of the atoms should be consistent with that in the simulation model file. The coordinates are in units of Ångstrom. This file is also generated by the hiphive package.