# Phono3py API#

How to use phono3py API is described below along with snippets that work with
`AlN-LDA`

example.

## Crystal structure#

Crystal structures in phono3py are usually `PhonopyAtoms`

class instances. When
we want to obtain the `PhonopyAtoms`

class instances from a file written in a
force-calculator format, the `read_crystal_structure`

function in phonopy may be
used, e.g., in the AlN-LDA example,

```
In [1]: from phonopy.interface.calculator import read_crystal_structure
In [2]: unitcell, _ = read_crystal_structure("POSCAR-unitcell", interface_mode='vasp')
In [3]: print(unitcell)
lattice:
- [ 3.110999999491908, 0.000000000000000, 0.000000000000000 ] # a
- [ -1.555499999745954, 2.694205030733368, 0.000000000000000 ] # b
- [ 0.000000000000000, 0.000000000000000, 4.978000000000000 ] # c
points:
- symbol: Al # 1
coordinates: [ 0.333333333333333, 0.666666666666667, 0.000948820000000 ]
mass: 26.981539
- symbol: Al # 2
coordinates: [ 0.666666666666667, 0.333333333333333, 0.500948820000000 ]
mass: 26.981539
- symbol: N # 3
coordinates: [ 0.333333333333333, 0.666666666666667, 0.619051180000000 ]
mass: 14.006700
- symbol: N # 4
coordinates: [ 0.666666666666667, 0.333333333333333, 0.119051180000000 ]
mass: 14.006700
```

Otherwise, it is directly created from
`PhonopyAtoms`

class.

```
In [1]: from phonopy.structure.atoms import PhonopyAtoms
In [2]: a = 3.11
In [3]: a = 3.111
In [4]: c = 4.978
In [5]: lattice = [[a, 0, 0], [-a / 2, a * np.sqrt(3) / 2, 0], [0, 0, c]]
In [6]: x = 1. / 3
In [7]: points = [[x, 2 * x, 0], [x * 2, x, 0.5], [x, 2 * x, 0.1181], [2 * x, x, 0.6181]]
In [8]: symbols = ['Al', 'Al', 'N', 'N']
In [9]: unitcell = PhonopyAtoms(cell=lattice, scaled_positions=points, symbols=symbols)
```

`Phono3py`

class#

To operate phono3py from python, there is phono3py API. The main class is
`Phono3py`

that is imported by

```
from phono3py import Phono3py
```

As written in Workflow, phono3py workflow is roughly divided into three
steps. `Phono3py`

class is used at each step. The minimum set of inputs to
instantiate `Phono3py`

class are `unitcell`

(1st argument), `supercell_matrix`

(2nd argument), and `primitive_matrix`

, which are used as

```
ph3 = Phono3py(unitcell, supercell_matrix=[3, 3, 2], primitive_matrix='auto')
```

The `unitcell`

is an instance of the
`PhonopyAtoms`

class.
`supercell_matrix`

and `primitive_matrix`

are the transformation matrices to
generate supercell and primitive cell from `unitcell`

(see
definitions).
This step is similar to the
instantiation of `Phonopy`

class.

There are many parameters that can be given to `Phono3py`

class. The details are
written in the docstring, which is shown by

```
help(Phono3py)
```

## Displacement dataset generation#

The step (1) in Workflow generates sets of displacements in supercell.
Supercells with the displacements are used as input crystal structure models of
force calculator that is, e.g., first-principles calculation code. Using the
force calculator, sets of forces of the supercells with the displacements are
obtained out of phono3py environment, i.e., phono3py only provides crystal
structures. Here we call the sets of the displacements as
`displacement dataset`

, the sets of the supercell forces as `force sets`

, and
the pair of `displacement dataset`

and `force sets`

as simply `dataset`

for
computing second and third force constants.

After instantiating `Phono3py`

with `unitcell`

, displacements are generated as
follows:

```
In [4]: ph3 = Phono3py(unitcell, supercell_matrix=[3, 3, 2], primitive_matrix='auto')
In [5]: ph3.generate_displacements()
In [6]: len(ph3.supercells_with_displacements)
Out[6]: 1254
In [7]: type(ph3.supercells_with_displacements[0])
Out[7]: phonopy.structure.atoms.PhonopyAtoms
```

By this, 1254 supercells with displacements were generated.

The generated displacement dataset is used in the force constants calculation.
Therefore, it is recommended to save it into a file if the python process having
the `Phono3py`

class instance is expected to be terminated. The displacement
dataset and crystal structure information are saved to a file by
`Phono3py.save()`

method:

```
In [8]: ph3.save("phono3py_disp.yaml")
```

## Supercell force calculation#

Forces of the generated supercells with displacements are calculated by some external force calculator such as first-principles calculation code.

Calculated supercell forces will be stored in a `Phono3py`

class instance
through `Phono3py.forces`

attribute by setting an array_like variable with the
shape of `(num_supercells, num_atoms_in_supercell, 3)`

. In the above example,
the array shape is `(1254, 72, 3)`

.

If the calculated force sets are stored in the
FORCES_FC3 file, the numpy array of `forces`

is
obtained by

```
forces = np.loadtxt("FORCES_FC3").reshape(-1, num_atoms_in_supercell, 3)
assert len(forces) == num_supercells
```

## Force constants calculation#

The pair of the displacement dataset and force sets is required to calculate
force constants. phono3py.load is the convenient function to load
these data from files and to set up them in the `Phono3py`

class instance.
However, in the case when only displacement dataset is expected, the low-level
phono3py-yaml parser in `Phono3pyYaml`

is useful.

```
In [1]: from phono3py.interface.phono3py_yaml import Phono3pyYaml
In [2]: ph3yml = Phono3pyYaml()
In [3]: ph3yml.read("phono3py_disp.yaml")
In [4]: disp_dataset = ph3yml.dataset
```

With this `ph3yml`

, how to compute force constants is explained. In the
following, it is assumed that we have `FORCES_FC3`

in the current directory. The
displacement dataset and force sets are set to the `Phono3py`

class instance as
follows:

```
In [5]: unitcell = ph3yml.unitcell
In [6]: from phono3py import Phono3py
In [7]: ph3 = Phono3py(unitcell, supercell_matrix=ph3yml.supercell_matrix, primitive_matrix=ph3yml.primitive_matrix)
In [8]: import numpy as np
In [9]: forces = np.loadtxt("FORCES_FC3").reshape(-1, len(ph3.supercell), 3)
In [10]: ph3.dataset = disp_dataset
In [11]: ph3.forces = forces
```

Now it is ready to compute force constants.

```
In [12]: ph3.produce_fc3()
```

## Non-analytical term correction parameters#

Users collects Born effective charges and dielectric constant tensor from experiments or calculations, and they are used as parameters for non-analytical term correction (NAC). These parameters are stored as a python dict:

```
'born': ndarray
Born effective charges
shape=(num_atoms_in_primitive_cell, 3, 3), dtype='double', order='C'
'factor': float
Unit conversion factor
'dielectric': ndarray
Dielectric constant tensor
shape=(3, 3), dtype='double', order='C'
```

Be careful that, in the case of using phono3py API, Born effective charges of
all atoms in the primitive cell are necessary, whereas in the
`BORN`

file,
only Born effective charges of symmetrically independent atoms are written. Some
more information about NAC parameters is found
here.

This NAC parameters are set to the `Phono3py`

class instance via the
`Phono3py.nac_params`

attribute. The NAC parameters may be read from `BORN`

file

```
In [13]: from phonopy.file_IO import parse_BORN
In [14]: nac_params = parse_BORN(ph3.primitive, filename="BORN")
In [15]: ph3.nac_params = nan_params
```

where the `parse_BORN`

function requires the corresponding primitive cell of the
`PhonopyAtoms`

class.

## Regular grid for **q**-point sampling#

Phonons are normally sampled on a \(\Gamma\) centre regular grid in reciprocal space. There are three ways to specify the regular grid.

Three integer values

One value

3x3 integer matrix

One of these is set to `mesh_numbers`

attribute of the `Phono3py`

class. For
example,

```
ph3.mesh_numbers = [10, 10, 10]
ph3.mesh_numbers = 50
ph3.mesh_numbers = [[-10, 10, 10], [10, -10, 10], [10, 10, -10]]
```

### Three integer values#

Three integer values (\(n_1\), \(n_2\), \(n_3\)) are specified.

The conventional regular grid is defined by the three integer values. The \(\mathbf{q}\)-points are given by linear combination of reciprocal basis vectors divided by the respective integers, which is given as

where \(m_i \in \{ 0, 1, \ldots, n_i - 1 \}\).

### One value#

A distance like value \(l\) is specified. Using this value, a regular grid is generated. As default, the three integer values of the conventional regular grid are defined by the following calculation.

Experimentally, use of a generalized regular grid is supported. By specifying
`use_grg = True`

at the `Phono3py`

class instantiation, the generalized regular
grid is generated using the value \(l\). First, the conventional unit cell of the
primitive cell is searched. Second, Eq. (2) is applied to the
reciprocal basis vectors of the conventional unit cell. Then,
\(\mathbf{q}\)-points are sampled following Eq. (1) with
respect to the reciprocal basis vectors of the conventional unit cell. The
parallelepiped defined by \(\mathbf{b}^*_i\) of the conventinal unit cell can be
smaller than that of the primitive cell. In this case, the \(\mathbf{q}\)-points
are sampled to fill the latter parallelepiped.

### 3x3 integer matrix (experimental)#

This is used to define the generalized regular grid explicitly. The generalized regular grid is designed to be automatically generated by one given value considering symmetry. However a 3x3 integer matrix (array) is accepted if this matrix follows the symmetry properly.

## Phonon-phonon interaction calculation#

Three phonon interaction strength is calculated after defining the regular grid.
When we have `phono3py_disp.yaml`

and `FORCES_FC3`

used above (maybe `BORN`

,
too), it is convenient to get the `Phono3py`

class instance with settings.

```
In [1]: import phono3py
In [2]: ph3 = phono3py.load("phono3py_disp.yaml", log_level=1)
NAC params were read from "BORN".
Displacement dataset for fc3 was read from "phono3py_disp.yaml".
Sets of supercell forces were read from "FORCES_FC3".
Computing fc3[ 1, x, x ] using numpy.linalg.pinv with displacements:
[ 0.0300 0.0000 0.0000]
[ 0.0000 0.0000 0.0300]
[ 0.0000 0.0000 -0.0300]
Computing fc3[ 37, x, x ] using numpy.linalg.pinv with displacements:
[ 0.0300 0.0000 0.0000]
[ 0.0000 0.0000 0.0300]
[ 0.0000 0.0000 -0.0300]
Expanding fc3.
fc3 was symmetrized.
Max drift of fc3: 0.000000 (yxx) 0.000000 (xyx) 0.000000 (xxy)
Displacement dataset for fc2 was read from "phono3py_disp.yaml".
Sets of supercell forces were read from "FORCES_FC3".
fc2 was symmetrized.
Max drift of fc2: 0.000000 (yy) 0.000000 (yy)
```

Three phonon interaction calculation becomes ready to run by the following lines:

```
In [3]: ph3.mesh_numbers = 30
In [4]: ph3.init_phph_interaction()
```

The three phonon interaction calculation is implemented in
`phono3py.phonon3.interaction.Interaction`

class. By phono3py’s design choice,
\(\mathbf{q}_1\) of the three phonons \((\mathbf{q}_1, \mathbf{q}_2, \mathbf{q}_3)\)
in this class is given, and calculation iterates over different \(\mathbf{q}_2\).
\(\mathbf{q}_3\) is uniquely determined from \(\mathbf{q}_1\) and \(\mathbf{q}_2\).
Symmetry is employed to avoid calculating symmetrically redundant \(\mathbf{q}_2\)
at the fixed \(\mathbf{q}_1\). One `Interaction`

class instance stores the data
only for one fixed \(\mathbf{q}_1\).

Three phonon interaction strength requires large memory space. Therefore, in lattice thermal conductivity calculation, it is calculated on-demand, and then abandoned after use of the data. Three phonon interaction strength calculation is the most computationally demanding part for usual applications. Detailed techniques are used to avoid elements of the Three phonon interaction strength if it is specified to do so.

Three phonon interaction strength calculation is the engine of more practical (or closer to macroscopic physical properties) calculation such as lattice thermal conductivity calculation. To minimize computational resources required by each purpose of use, the outer function that calls this function determines proper computational configuration of this function.

## Lattice thermal conductivity calculation#

Once three phonon interaction calculation is prepared by
`ph3.init_phph_interaction()`

, it is ready to run lattice thermal conductivity
calculation. The many parameters are explained in the docstring, but even with
the default parameters, the lattice thermal conductivity calculation under the
mode relaxation time approximation is performed as follows:

```
In [5]: ph3.run_thermal_conductivity()
```

## Use of different supercell dimensions for 2nd and 3rd order FCs#

Phono3py supports different supercell dimensions for second and third order
force constants (fc2 and fc3). The default setting is using the same dimensions.
Although fc3 requires much more number of supercells than fc2, in our
experience, fc3 have shorter interaction range in direct space. Therefore we
tend to expect to use smaller supercell dimension for fc3 than that for fc2. To
achieve this, `phonon_supercell_matrix`

parameter exists to specify fc2
supercell dimension independently.

Using `POSCAR-unitcell`

in the AlN-LDA example,

```
In [1]: from phonopy.interface.calculator import read_crystal_structure
In [2]: unitcell, _ = read_crystal_structure("POSCAR-unitcell", interface_mode='vasp')
In [3]: from phono3py import Phono3py
In [4]: ph3 = Phono3py(unitcell, supercell_matrix=[3, 3, 2], primitive_matrix='auto', phonon_supercell_matrix=[5, 5, 3])
In [5]: ph3.save("phono3py_disp.yaml")
In [6]: ph3.generate_displacements()
```

`Phono3py.generate_fc2_displacements()`

is the method to generate displacement
dataset for fc2. This is normally unnecessary to be called because this is
called when `Phono3py.generate_displacements()`

is called.
`Phono3py.generate_fc2_displacements()`

may be called explicitly if non-default
control of displacement pattern is expected. See their docstrings for the
details.

When `phonon_supercell_matrix`

is specified, the following attributes are
usable:

```
Phono3py.phonon_supercell_matrix
Phono3py.phonon_dataset
Phono3py.phonon_forces
Phono3py.phonon_supercells_with_displacements
```

The meanings of them are found in their docstrings though they may be guessed easily.

`phono3py.load`

#

The purpose of `phono3py.load`

is to create a `Phono3py`

class instance with
basic parameters loaded from the phono3py-yaml-type file and also to try setting
up force constants and non-analytical term correction automatically from
phono3py files in the current directory.

In AlN-LDA example, the unit cell structure, supercell matrix, and primitive
matrix were recorded in the `phono3py_disp_dimfc2.yaml`

file. This is easily
read a helper function of `phono3py.load`

. Using ipython (or jupyter-notebook):

```
In [1]: import phono3py
In [2]: ph3 = phono3py.load("phono3py_disp_dimfc2.yaml", produce_fc=False)
In [3]: type(ph3)
Out[3]: phono3py.api_phono3py.Phono3py
In [4]: print(ph3.unitcell)
lattice:
- [ 3.110999999491908, 0.000000000000000, 0.000000000000000 ] # a
- [ -1.555499999745954, 2.694205030733368, 0.000000000000000 ] # b
- [ 0.000000000000000, 0.000000000000000, 4.978000000000000 ] # c
points:
- symbol: Al # 1
coordinates: [ 0.333333333333333, 0.666666666666667, 0.000948820000000 ]
mass: 26.981539
- symbol: Al # 2
coordinates: [ 0.666666666666667, 0.333333333333333, 0.500948820000000 ]
mass: 26.981539
- symbol: N # 3
coordinates: [ 0.333333333333333, 0.666666666666667, 0.619051180000000 ]
mass: 14.006700
- symbol: N # 4
coordinates: [ 0.666666666666667, 0.333333333333333, 0.119051180000000 ]
mass: 14.006700
In [5]: ph3.supercell_matrix
Out[5]:
array([[3, 0, 0],
[0, 3, 0],
[0, 0, 2]])
In [6]: ph3.nac_params
Out[6]:
{'born': array([[[ 2.51264750e+00, 0.00000000e+00, 0.00000000e+00],
[ 0.00000000e+00, 2.51264750e+00, 0.00000000e+00],
[ 0.00000000e+00, 0.00000000e+00, 2.67353000e+00]],
[[ 2.51264750e+00, -2.22044605e-16, 0.00000000e+00],
[ 0.00000000e+00, 2.51264750e+00, 0.00000000e+00],
[ 0.00000000e+00, 0.00000000e+00, 2.67353000e+00]],
[[-2.51264750e+00, 0.00000000e+00, 0.00000000e+00],
[ 0.00000000e+00, -2.51264750e+00, 0.00000000e+00],
[ 0.00000000e+00, 0.00000000e+00, -2.67353000e+00]],
[[-2.51264750e+00, 2.22044605e-16, 0.00000000e+00],
[ 0.00000000e+00, -2.51264750e+00, 0.00000000e+00],
[ 0.00000000e+00, 0.00000000e+00, -2.67353000e+00]]]),
'factor': 14.39965172592227,
'dielectric': array([[4.435009, 0. , 0. ],
[0. , 4.435009, 0. ],
[0. , 0. , 4.653269]])}
```