mirror of https://github.com/phonopy/phonopy.git
769 lines
23 KiB
Markdown
769 lines
23 KiB
Markdown
(phonopy_module)=
|
|
|
|
# Phonopy API for Python
|
|
|
|
**This is under development. Configurations may alter.** Requests or suggestions
|
|
are very welcome.
|
|
|
|
## Import modules
|
|
|
|
After setting the phonopy python path, the phonopy module is imported by:
|
|
|
|
```python
|
|
from phonopy import Phonopy
|
|
```
|
|
|
|
Crystal structure is defined by the `PhonopyAtoms` class. The `PhonopyAtoms`
|
|
module is imported by:
|
|
|
|
```python
|
|
from phonopy.structure.atoms import PhonopyAtoms
|
|
```
|
|
|
|
The instance of `PhonopyAtoms` can be made by reading a crystal structure in a
|
|
variety of calculator formats found at {ref}`calculator_interfaces`.
|
|
|
|
```python
|
|
from phonopy.interface.calculator import read_crystal_structure
|
|
unitcell, _ = read_crystal_structure("POSCAR-unitcell", interface_mode='vasp')
|
|
```
|
|
|
|
For VASP format, the keyword argument of `interface_mode` can be omitted. For
|
|
QE,
|
|
|
|
```python
|
|
unitcell, optional_structure_info = read_crystal_structure("NaCl.in", interface_mode='qe')
|
|
```
|
|
|
|
Note that `read_crystal_structure` returns a tuple and the first element is the
|
|
`PhonopyAtoms` instance.
|
|
## Work flow
|
|
|
|
The work flow is schematically shown in {ref}`workflow`.
|
|
|
|
### Pre-process
|
|
|
|
The first step is to create a `Phonopy` object with at least two arguments, a
|
|
unit cell (`PhonopyAtoms` object, see {ref}`phonopy_Atoms`) and a supercell
|
|
matrix (3x3 array, see {ref}`variable_supercell_matrix`). In the following
|
|
example, a {math}`2\times 2\times 2` supercell is created. The displacements to
|
|
be introduced to the supercell are internally generated by the
|
|
`generate_displacements()` method with the `distance` keyword argument. The
|
|
supercells with displacements are obtained by
|
|
`get_supercells_with_displacements()` method as a list of `PhonopyAtoms`
|
|
objects.
|
|
|
|
```python
|
|
import numpy as np
|
|
from phonopy import Phonopy
|
|
from phonopy.structure.atoms import PhonopyAtoms
|
|
|
|
a = 5.404
|
|
unitcell = PhonopyAtoms(symbols=['Si'] * 8,
|
|
cell=(np.eye(3) * a),
|
|
scaled_positions=[[0, 0, 0],
|
|
[0, 0.5, 0.5],
|
|
[0.5, 0, 0.5],
|
|
[0.5, 0.5, 0],
|
|
[0.25, 0.25, 0.25],
|
|
[0.25, 0.75, 0.75],
|
|
[0.75, 0.25, 0.75],
|
|
[0.75, 0.75, 0.25]])
|
|
phonon = Phonopy(unitcell,
|
|
supercell_matrix=[[2, 0, 0], [0, 2, 0], [0, 0, 2]],
|
|
primitive_matrix=[[0, 0.5, 0.5],
|
|
[0.5, 0, 0.5],
|
|
[0.5, 0.5, 0]])
|
|
phonon.generate_displacements(distance=0.03)
|
|
supercells = phonon.supercells_with_displacements
|
|
```
|
|
|
|
In this example, the displacement distance is set to 0.03 (in Angstrom if the
|
|
crystal structure uses the Angstrom unit and the default value is 0.01.) The
|
|
supercells with displacements are given as a list of `PhonopyAtoms`. See
|
|
{ref}`phonopy_read_write_structure` to write
|
|
those into files in a crystal structure format.
|
|
|
|
The frequency unit conversion factor to THz has to be set by using the `factor`
|
|
keyword in `Phonopy` class. The factors are `VaspToTHz` for VASP, `Wien2kToTHz`
|
|
for Wien2k, `AbinitToTHz` for Abinit, `PwscfToTHz` for Pwscf, `ElkToTHz` for
|
|
Elk, `SiestaToTHz` for Siesta, `CrystalToTHz` for CRYSTAL, `FleurToTHz` for
|
|
Fleur, `VaspToTHz`, and `DftbpToTHz` for DFTB+ is the default value. For
|
|
example:
|
|
|
|
```python
|
|
from phonopy.units import AbinitToTHz
|
|
|
|
phonon = Phonopy(unitcell,
|
|
supercell_matrix=[[2, 0, 0], [0, 2, 0], [0, 0, 2]],
|
|
primitive_matrix=[[0, 0.5, 0.5],
|
|
[0.5, 0, 0.5],
|
|
[0.5, 0.5, 0]],
|
|
factor=AbinitToTHz)
|
|
```
|
|
|
|
Some more information on physical unit conversion is found at
|
|
{ref}`frequency_conversion_factor_tag`, {ref}`physical_unit_conversion`, and
|
|
{ref}`calculator_interfaces`.
|
|
|
|
### Post process
|
|
|
|
Forces on atoms are supposed to be obtained by running force calculator (e.g.
|
|
VASP) with each supercell with a displacement. Then the forces in the
|
|
calculation outputs have to be collected by users. However output parsers for
|
|
selected calculators are found under `phonopy.interface`, which may be useful.
|
|
The forces have to be stored in a specific structure: a numpy array (or nested
|
|
list) as follows:
|
|
|
|
```python
|
|
[ [ [ f_1x, f_1y, f_1z ], [ f_2x, f_2y, f_2z ], ... ], # first supercell
|
|
[ [ f_1x, f_1y, f_1z ], [ f_2x, f_2y, f_2z ], ... ], # second supercell
|
|
... ]
|
|
```
|
|
|
|
This array (`sets_of_forces`) is set to the `Phonopy` object by:
|
|
|
|
```python
|
|
phonon.forces = sets_of_forces
|
|
```
|
|
|
|
This is the case when the set of atomic displacements is generated internally.
|
|
The information of displacements is already stored in the `Phonopy` object. But
|
|
if you want to input the forces together with the corresponding custom set of
|
|
displacements, `displacement_dataset` has to be prepared as a python dictionary
|
|
as follows:
|
|
|
|
```python
|
|
displacement_dataset =
|
|
{'natom': number_of_atoms_in_supercell,
|
|
'first_atoms': [
|
|
{'number': atom index of displaced atom (starting with 0),
|
|
'displacement': displacement in Cartesian coordinates,
|
|
'forces': forces on atoms in supercell},
|
|
{...}, ...]}
|
|
```
|
|
|
|
This is set to the `Phonopy` object by:
|
|
|
|
```python
|
|
phonon.dataset = displacement_dataset
|
|
```
|
|
|
|
From the set of displacements and forces, force constants internally with
|
|
calculated supercell sets of forces by
|
|
|
|
```python
|
|
phonon.produce_force_constants()
|
|
```
|
|
|
|
If you have force constants and don't need to create force constants from forces
|
|
and displacements, simply set your force constants by
|
|
|
|
```python
|
|
phonon.force_constants = force_constants
|
|
```
|
|
|
|
The force constants matrix is given in 4 dimensional array (better to be a numpy
|
|
array of `dtype='double', order='C'`). The shape of force constants matrix is
|
|
`(N, N, 3, 3)` where `N` is the number of atoms in the supercell and 3 gives
|
|
Cartesian axes. The compact force constants matrix with `(Np, N, 3, 3)` where
|
|
`Np` is the number of atoms in the primitive cell is also supported. See the
|
|
details at {ref}`file_force_constants`.
|
|
|
|
## Phonon calculation
|
|
|
|
(phonopy_save_parameters)=
|
|
|
|
### Save parameters (`phonopy.save`)
|
|
|
|
Basic information and parameters needed for phonon calculation are saved into a
|
|
file by `phonopy.save`.
|
|
|
|
```python
|
|
phonon.save()
|
|
```
|
|
|
|
The default file name is `phonopy_params.yaml`. Force sets, displacements, Born
|
|
effective charges, and dielectric constant are written in the default behaviour.
|
|
If force constants are needed to be written in the yaml file, the argument
|
|
`settings` is set as follows:
|
|
|
|
```python
|
|
phonon.save(settings={'force_constants': True})
|
|
```
|
|
|
|
|
|
### Band structure
|
|
|
|
Set band paths (`run_band_structure()`) and get the results
|
|
(`get_band_structure_dict()`).
|
|
|
|
A dictionary with `qpoints`, `distances`, `frequencies`, `eigenvectors`,
|
|
`group_velocities` is returned by `get_band_structure_dict()`. Eigenvectors can
|
|
be obtained when `with_eigenvectors=True` at `run_band_structure()`. See the
|
|
details at docstring of `Phonopy.get_band_structure_dict`. Phonon frequency is
|
|
sqrt(eigenvalue). A negative eigenvalue has to correspond to the imaginary
|
|
frequency, but for the plotting, it is set as the negative value in the above
|
|
example. In addition, you need to multiply by your unit conversion factor. In
|
|
the case of VASP to transform to THz, the factor is 15.633302.
|
|
|
|
In `example/NaCl`, the phonopy is executed from python script, e.g.,
|
|
|
|
```python
|
|
import phonopy
|
|
from phonopy.phonon.band_structure import get_band_qpoints_and_path_connections
|
|
|
|
path = [[[0, 0, 0], [0.5, 0, 0.5], [0.625, 0.25, 0.625]],
|
|
[[0.375, 0.375, 0.75], [0, 0, 0], [0.5, 0.5, 0.5], [0.5, 0.25, 0.75]]]
|
|
labels = ["$\\Gamma$", "X", "U", "K", "$\\Gamma$", "L", "W"]
|
|
qpoints, connections = get_band_qpoints_and_path_connections(path, npoints=51)
|
|
phonon = phonopy.load("phonopy_disp.yaml")
|
|
phonon.run_band_structure(qpoints, path_connections=connections, labels=labels)
|
|
phonon.plot_band_structure().show()
|
|
|
|
# To plot DOS next to band structure
|
|
phonon.run_mesh([20, 20, 20])
|
|
phonon.run_total_dos()
|
|
phonon.plot_band_structure_and_dos().show()
|
|
|
|
# To plot PDOS next to band structure
|
|
phonon.run_mesh([20, 20, 20], with_eigenvectors=True, is_mesh_symmetry=False)
|
|
phonon.run_projected_dos()
|
|
phonon.plot_band_structure_and_dos(pdos_indices=[[0], [1]]).show()
|
|
```
|
|
|
|
`path_connections` and `labels` are unnecessary to set unless nice looking
|
|
plotting is needed. To obtain eigenvectors, it is necessary to inform to store
|
|
eigenvectors by:
|
|
|
|
```python
|
|
phonon.run_band_structure(bands, with_eigenvectors=True)
|
|
```
|
|
|
|
To obtain group velocities:
|
|
|
|
```python
|
|
phonon.run_band_structure(bands, with_group_velocities=True)
|
|
```
|
|
|
|
Automatic selection of band paths using
|
|
[SeeK-path](https://seekpath.readthedocs.io/en/latest/) is invoked by
|
|
|
|
```python
|
|
phonon.auto_band_structure()
|
|
```
|
|
|
|
and to plot
|
|
|
|
```python
|
|
phonon.auto_band_structure(plot=True).show()
|
|
```
|
|
|
|
To use this method, `seekpath` python module is needed.
|
|
|
|
### Mesh sampling
|
|
|
|
Set sampling mesh (`set_mesh`) in reciprocal space. The irreducible _q_-points
|
|
and corresponding _q_-point weights, eigenvalues, and eigenvectors are obtained
|
|
by `get_mesh_dict()`. `mesh` gives the sampling mesh with Monkhorst-Pack scheme.
|
|
The keyword `shift` gives the fractional mesh shift with respect to the
|
|
neighboring grid points.
|
|
|
|
```python
|
|
mesh = [20, 20, 20]
|
|
phonon.run_mesh(mesh)
|
|
mesh_dict = phonon.get_mesh_dict()
|
|
qpoints = mesh_dict['qpoints']
|
|
weights = mesh_dict['weights']
|
|
frequencies = mesh_dict['frequencies']
|
|
eigenvectors = mesh_dict['eigenvectors']
|
|
group_velocities = mesh_dict['group_velocities']
|
|
```
|
|
|
|
To obtain eigenvectors, it is necessary to inform to store eigenvectors by:
|
|
|
|
```python
|
|
phonon.run_mesh([20, 20, 20], with_eigenvectors=True)
|
|
```
|
|
|
|
and for group velocities:
|
|
|
|
```python
|
|
phonon.run_mesh([20, 20, 20], with_group_velocities=True)
|
|
```
|
|
|
|
The first argument of `run_mesh()` can be a float value, which is a length
|
|
measure as explained at {ref}`mesh_tag`, for example:
|
|
|
|
```python
|
|
phonon.run_mesh(100.0)
|
|
```
|
|
|
|
### DOS and PDOS
|
|
|
|
Before starting mesh sampling has to be finished. Then set parameters
|
|
(`run_total_dos()` or `run_projected_dos()`) and write the results into files
|
|
(`write_total_dos()` and `write_projected_dos()`). In the case of PDOS, the
|
|
eigenvectors have to be calculated in the mesh sampling. To get the results
|
|
`get_total_dos_dict()` and `get_projected_dos_dict()` can be used.
|
|
|
|
To plot total DOS,
|
|
|
|
```python
|
|
phonon.run_mesh([20, 20, 20]) phonon.run_total_dos()
|
|
phonon.plot_total_dos().show()
|
|
```
|
|
|
|
and projected DOS
|
|
|
|
```python
|
|
phonon.run_mesh([20, 20, 20], with_eigenvectors=True, is_mesh_symmetry=False)
|
|
phonon.run_projected_dos()
|
|
phonon.plot_projected_dos().show()
|
|
```
|
|
|
|
Convenient shortcuts exist as follows:
|
|
|
|
```python
|
|
phonon.auto_total_dos(plot=True).show()
|
|
```
|
|
|
|
and
|
|
|
|
```python
|
|
phonon.auto_projected_dos(plot=True).show()
|
|
```
|
|
|
|
### Thermal properties
|
|
|
|
Before starting the thermal property calculation, the mesh sampling calculation
|
|
has to be done in the **THz unit**. The unit conversion factor for phonon
|
|
frequency is set in the pre-process of Phonopy with the `factor` keyword.
|
|
Calculation range of temperature is set by the parameters
|
|
`run_thermal_properties`. Helmholtz free energy, entropy, heat capacity at
|
|
constant volume at temperatures are obtained by `get_thermal_properties_dict`,
|
|
where the results are given as a dictionary of temperatures, Helmholtz free
|
|
energy, entropy, and heat capacity with keys `temperatures`, `free_energy`,
|
|
`entropy`, and `heat_capacity`, respectively.
|
|
|
|
```python
|
|
phonon.run_mesh([20, 20, 20])
|
|
phonon.run_thermal_properties(t_step=10,
|
|
t_max=1000,
|
|
t_min=0)
|
|
tp_dict = phonon.get_thermal_properties_dict()
|
|
temperatures = tp_dict['temperatures']
|
|
free_energy = tp_dict['free_energy']
|
|
entropy = tp_dict['entropy']
|
|
heat_capacity = tp_dict['heat_capacity']
|
|
|
|
for t, F, S, cv in zip(temperatures, free_energy, entropy, heat_capacity):
|
|
print(("%12.3f " + "%15.7f" * 3) % ( t, F, S, cv ))
|
|
|
|
phonon.plot_thermal_properties().show()
|
|
```
|
|
|
|
### Non-analytical term correction
|
|
|
|
To apply non-analytical term correction, Born effective charge tensors for all
|
|
atoms in **primitive** cell, dielectric constant tensor, and the unit conversion
|
|
factor have to be correctly set. The tensors are given in Cartesian coordinates.
|
|
|
|
```python
|
|
born = [[[1.08878299, 0, 0],
|
|
[0, 1.08878299, 0],
|
|
[0, 0, 1.08878299]],
|
|
[[-1.08878299, 0, 0],
|
|
[0, -1.08878299, 0],
|
|
[0, 0, -1.08878299]]]
|
|
epsilon = [[2.56544559, 0, 0],
|
|
[0, 2.56544559, 0],
|
|
[0, 0, 2.56544559]]
|
|
factors = 14.400
|
|
phonon.nac_params = {'born': born,
|
|
'factor': factors,
|
|
'dielectric': epsilon}
|
|
```
|
|
|
|
## Data structure
|
|
|
|
### Eigenvectors
|
|
|
|
Eigenvectors are given as the column vectors. Internally phonopy uses
|
|
`numpy.linalg.eigh` and `eigh` is a wrapper of LAPACK. So eigenvectors follow
|
|
the convention of LAPACK, which can be shown at
|
|
http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eigh.html
|
|
|
|
Eigenvectors corresponding to phonopy yaml output are obtained as follows.
|
|
|
|
#### Band structure
|
|
|
|
```python
|
|
if eigvecs is not None:
|
|
for eigvecs_on_path in eigvecs:
|
|
for eigvecs_at_q in eigvecs_on_path:
|
|
for vec in eigvecs_at_q.T:
|
|
print(vec)
|
|
```
|
|
|
|
#### Mesh sampling
|
|
|
|
```python
|
|
if eigvecs is not None:
|
|
for eigvecs_at_q in eigvecs:
|
|
for vec in eigvecs_at_q.T:
|
|
print(vec)
|
|
```
|
|
|
|
|
|
(phonopy_Atoms)=
|
|
## `PhonopyAtoms` class
|
|
|
|
### Initialization
|
|
|
|
The usable keywords in the initialization are:
|
|
|
|
```python
|
|
cell=None,
|
|
scaled_positions=None,
|
|
positions=None,
|
|
numbers=None,
|
|
symbols=None,
|
|
masses=None,
|
|
magnetic_moments=None,
|
|
```
|
|
|
|
At least three arguments have to be given at the initialization, which are
|
|
|
|
- `cell`
|
|
- `positions` or `scaled_positions`
|
|
- `symbols` or `numbers`
|
|
|
|
(phonopy_Atoms_variables)=
|
|
|
|
### Variables
|
|
|
|
The following variables are implemented in the `PhonopyAtoms` class in
|
|
`phonopy/structure/atoms.py`.
|
|
|
|
(phonopy_Atoms_cell)=
|
|
|
|
#### `cell`
|
|
|
|
Basis vectors are given in the matrix form in Cartesian coordinates.
|
|
|
|
```python
|
|
[ [ a_x, a_y, a_z ], [ b_x, b_y, b_z ], [ c_x, c_y, c_z ] ]
|
|
```
|
|
|
|
#### `scaled_positions`
|
|
|
|
Atomic positions in fractional coordinates.
|
|
|
|
```python
|
|
[ [ x1_a, x1_b, x1_c ], [ x2_a, x2_b, x2_c ], [ x3_a, x3_b, x3_c ], ... ]
|
|
```
|
|
|
|
#### `positions`
|
|
|
|
Cartesian positions of atoms.
|
|
|
|
```python
|
|
positions = np.dot(scaled_positions, cell)
|
|
```
|
|
|
|
where `np` means the numpy module (`import numpy as np`).
|
|
|
|
#### `symbols`
|
|
|
|
Chemical symbols, e.g.,
|
|
|
|
```python
|
|
['Zn', 'Zn', 'O', 'O']
|
|
```
|
|
|
|
for the ZnO unit cell.
|
|
|
|
#### `numbers`
|
|
|
|
Atomic numbers, e.g.,
|
|
|
|
```python
|
|
[30, 30, 8, 8]
|
|
```
|
|
|
|
for the ZnO unit cell.
|
|
|
|
#### `masses`
|
|
|
|
Atomic masses, e.g.,
|
|
|
|
```python
|
|
[65.38, 65.38, 15.9994, 15.9994]
|
|
```
|
|
|
|
for the ZnO unit cell.
|
|
|
|
### Attributes
|
|
|
|
```
|
|
cell
|
|
positions
|
|
scaled_positions
|
|
masses
|
|
magnetic_moments
|
|
symbols
|
|
numbers
|
|
volume
|
|
```
|
|
|
|
where `volume` is the getter only.
|
|
|
|
### Methods
|
|
|
|
`unitcell.get_number_of_atoms()` is equivalent to `len(unitcell)`. An instance
|
|
can be deep-copied by `unitcell.copy()`. Human-readable crystal structure in
|
|
Yaml format is shown by `print(unitcell)`. `unitcell.to_tuple` converts to
|
|
spglib crystal structure
|
|
(https://spglib.github.io/spglib/python-spglib.html#crystal-structure-cell).
|
|
|
|
## Definitions of variables
|
|
|
|
(variable_primitive_matrix)=
|
|
|
|
### Primitive matrix
|
|
|
|
Primitive matrix {math}`M_\mathrm{p}` is a tranformation matrix from lattice
|
|
vectors to those of a primitive cell if there exists the primitive cell in the
|
|
lattice vectors. Following a crystallography convention, the transformation is
|
|
given by
|
|
|
|
```{math}
|
|
( \mathbf{a}_\mathrm{p} \; \mathbf{b}_\mathrm{p} \; \mathbf{c}_\mathrm{p} ) = (
|
|
\mathbf{a}_\mathrm{u} \; \mathbf{b}_\mathrm{u} \; \mathbf{c}_\mathrm{u} )
|
|
M_\mathrm{p}
|
|
```
|
|
|
|
where {math}`\mathbf{a}_\mathrm{u}`, {math}`\mathbf{b}_\mathrm{u}`, and
|
|
{math}`\mathbf{c}_\mathrm{u}` are the column vectors of the original lattice
|
|
vectors, and {math}`\mathbf{a}_\mathrm{p}`, {math}`\mathbf{b}_\mathrm{p}`, and
|
|
{math}`\mathbf{c}_\mathrm{p}` are the column vectors of the primitive lattice
|
|
vectors. Be careful that the lattice vectors of the `PhonopyAtoms` class are the
|
|
row vectors ({ref}`phonopy_Atoms_cell`). Therefore the phonopy code, which
|
|
relies on the `PhonopyAtoms` class, is usually written such as
|
|
|
|
```python
|
|
primitive_lattice = np.dot(original_lattice.T, primitive_matrix).T,
|
|
```
|
|
|
|
or equivalently,
|
|
|
|
```python
|
|
primitive_lattice = np.dot(primitive_matrix.T, original_lattice)
|
|
```
|
|
|
|
(variable_supercell_matrix)=
|
|
|
|
### Supercell matrix
|
|
|
|
Supercell matrix {math}`M_\mathrm{s}` is a transformation matrix from lattice
|
|
vectors to those of a super cell. Following a crystallography convention, the
|
|
transformation is given by
|
|
|
|
```{math}
|
|
( \mathbf{a}_\mathrm{s} \; \mathbf{b}_\mathrm{s} \; \mathbf{c}_\mathrm{s} ) = (
|
|
\mathbf{a}_\mathrm{u} \; \mathbf{b}_\mathrm{u} \; \mathbf{c}_\mathrm{u} )
|
|
M_\mathrm{s}
|
|
```
|
|
|
|
where {math}`\mathbf{a}_\mathrm{u}`, {math}`\mathbf{b}_\mathrm{u}`, and
|
|
{math}`\mathbf{c}_\mathrm{u}` are the column vectors of the original lattice
|
|
vectors, and {math}`\mathbf{a}_\mathrm{s}`, {math}`\mathbf{b}_\mathrm{s}`, and
|
|
{math}`\mathbf{c}_\mathrm{s}` are the column vectors of the supercell lattice
|
|
vectors. Be careful that the lattice vectors of the `PhonopyAtoms` class are the
|
|
row vectors ({ref}`phonopy_Atoms_cell`). Therefore the phonopy code, which
|
|
relies on the `PhonopyAtoms` class, is usually written such as
|
|
|
|
```python
|
|
supercell_lattice = np.dot(original_lattice.T, supercell_matrix).T,
|
|
```
|
|
|
|
or equivalently,
|
|
|
|
```python
|
|
supercell_lattice = np.dot(supercell_matrix.T, original_lattice)
|
|
```
|
|
|
|
### Symmetry search tolerance
|
|
|
|
Symmetry search tolerance (often the name `symprec` is used in phonopy) is used
|
|
to determine symmetry operations of the crystal structures. The physical unit
|
|
follows that of input crystal structure.
|
|
|
|
(phonopy_load)=
|
|
## Load phonopy settings `phonopy.load`
|
|
|
|
`phonopy.load` is a convenient function that creates a `Phonopy` instance by
|
|
loading data from a `phonopy_xxx.yaml` file, which may include all the necessary
|
|
information to run phonopy. A typical usage is:
|
|
|
|
```python
|
|
import phonopy
|
|
phonon = phonopy.load("phonopy_params.yaml")
|
|
```
|
|
|
|
If `phonopy_params.yaml` contains a displacement-force dataset and you want to
|
|
avoid producing force constants, set `produce_fc=False`:
|
|
|
|
```python
|
|
phonon = phonopy.load("phonopy_params.yaml", produce_fc=False)
|
|
```
|
|
|
|
Alternatively, if you have either `phonopy.yaml` (or `phonopy_disp.yaml`) along
|
|
with a `FORCE_SETS` file, you can create a `Phonopy` object like this:
|
|
|
|
```python
|
|
phonon = phonopy.load("phonopy.yaml", force_sets_filename="FORCE_SETS")
|
|
```
|
|
|
|
In this case, the command reads the structure information from `phonopy.yaml`
|
|
and the displacement-force data from `FORCE_SETS` to create the `Phonopy`
|
|
instance.
|
|
|
|
If your current directory contains the following files:
|
|
|
|
```bash
|
|
% ls
|
|
BORN FORCE_SETS phonopy.yaml
|
|
```
|
|
|
|
then both the `BORN` and `FORCE_SETS` files will be read automatically by
|
|
|
|
```python
|
|
phonon = phonopy.load("phonopy.yaml")
|
|
```
|
|
|
|
For more details, see the function's docstring:
|
|
|
|
```python
|
|
In [1]: import phonopy
|
|
In [2]: help(phonopy.load)
|
|
```
|
|
|
|
(phonopy_read_write_structure)=
|
|
## Read and write crystal structures
|
|
|
|
There is a function to write the `PhonopyAtoms` instance into crystal structure
|
|
formats of different force calculators, `write_crystal_structure`. This works as
|
|
a partner of `read_crystal_structure`. Taking an example of QE interface, how to
|
|
use these functions is shown below.
|
|
|
|
```ipython
|
|
In [1]: from phonopy.interface.calculator import read_crystal_structure, write_crystal_structure
|
|
|
|
In [2]: !cat "NaCl.in"
|
|
&control
|
|
calculation = 'scf'
|
|
tprnfor = .true.
|
|
tstress = .true.
|
|
pseudo_dir = '/home/togo/espresso/pseudo/'
|
|
/
|
|
&system
|
|
ibrav = 0
|
|
nat = 8
|
|
ntyp = 2
|
|
ecutwfc = 70.0
|
|
/
|
|
&electrons
|
|
diagonalization = 'david'
|
|
conv_thr = 1.0d-9
|
|
/
|
|
ATOMIC_SPECIES
|
|
Na 22.98976928 Na.pbe-spn-kjpaw_psl.0.2.UPF
|
|
Cl 35.453 Cl.pbe-n-kjpaw_psl.0.1.UPF
|
|
ATOMIC_POSITIONS crystal
|
|
Na 0.0000000000000000 0.0000000000000000 0.0000000000000000
|
|
Na 0.0000000000000000 0.5000000000000000 0.5000000000000000
|
|
Na 0.5000000000000000 0.0000000000000000 0.5000000000000000
|
|
Na 0.5000000000000000 0.5000000000000000 0.0000000000000000
|
|
Cl 0.5000000000000000 0.5000000000000000 0.5000000000000000
|
|
Cl 0.5000000000000000 0.0000000000000000 0.0000000000000000
|
|
Cl 0.0000000000000000 0.5000000000000000 0.0000000000000000
|
|
Cl 0.0000000000000000 0.0000000000000000 0.5000000000000000
|
|
CELL_PARAMETERS angstrom
|
|
5.6903014761756712 0 0
|
|
0 5.6903014761756712 0
|
|
0 0 5.6903014761756712
|
|
K_POINTS automatic
|
|
8 8 8 1 1 1
|
|
|
|
In [3]: cell, optional_structure_info = read_crystal_structure("NaCl.in", interface_mode='qe')
|
|
|
|
In [4]: optional_structure_info
|
|
Out[4]:
|
|
('NaCl.in',
|
|
{'Na': 'Na.pbe-spn-kjpaw_psl.0.2.UPF', 'Cl': 'Cl.pbe-n-kjpaw_psl.0.1.UPF'})
|
|
|
|
In [5]: write_crystal_structure("NaCl-out.in", cell, interface_mode='qe', optional_structure_info=optional_structure_info)
|
|
|
|
In [6]: !cat "NaCl-out.in"
|
|
! ibrav = 0, nat = 8, ntyp = 2
|
|
CELL_PARAMETERS bohr
|
|
10.7531114272216008 0.0000000000000000 0.0000000000000000
|
|
0.0000000000000000 10.7531114272216008 0.0000000000000000
|
|
0.0000000000000000 0.0000000000000000 10.7531114272216008
|
|
ATOMIC_SPECIES
|
|
Na 22.98977 Na.pbe-spn-kjpaw_psl.0.2.UPF
|
|
Cl 35.45300 Cl.pbe-n-kjpaw_psl.0.1.UPF
|
|
ATOMIC_POSITIONS crystal
|
|
Na 0.0000000000000000 0.0000000000000000 0.0000000000000000
|
|
Na 0.0000000000000000 0.5000000000000000 0.5000000000000000
|
|
Na 0.5000000000000000 0.0000000000000000 0.5000000000000000
|
|
Na 0.5000000000000000 0.5000000000000000 0.0000000000000000
|
|
Cl 0.5000000000000000 0.5000000000000000 0.5000000000000000
|
|
Cl 0.5000000000000000 0.0000000000000000 0.0000000000000000
|
|
Cl 0.0000000000000000 0.5000000000000000 0.0000000000000000
|
|
Cl 0.0000000000000000 0.0000000000000000 0.5000000000000000
|
|
```
|
|
|
|
Depending on calculator interfaces, all the information can not be recovered
|
|
from the information obtained from `read_crystal_structure`. More details about
|
|
how `write_crystal_structure` works may need to read directly the
|
|
[code](https://github.com/phonopy/phonopy/blob/develop/phonopy/interface/calculator.py#L123).
|
|
|
|
## Getting parameters for non-analytical term correction
|
|
|
|
Parameters for non-analytical term correction may be made as follows. This
|
|
example assumes that the user knows what are the unit cell and primitive cell
|
|
and that the Born effective charge and dielectric constant were calculated using
|
|
VASP code by the unit cell.
|
|
|
|
```python
|
|
import io
|
|
import numpy as np
|
|
from phonopy.units import Hartree, Bohr
|
|
from phonopy.structure.symmetry import symmetrize_borns_and_epsilon
|
|
from phonopy.interface.vasp import VasprunxmlExpat
|
|
|
|
with io.open("vasprun.xml", "rb") as f:
|
|
vasprun = VasprunxmlExpat(f)
|
|
vasprun.parse():
|
|
epsilon = vasprun.epsilon
|
|
borns = vasprun.born
|
|
unitcell = vasprun.cell
|
|
|
|
borns_, epsilon_ = symmetrize_borns_and_epsilon(
|
|
borns,
|
|
epsilon,
|
|
unitcell,
|
|
primitive_matrix=[[0, 0.5, 0.5],
|
|
[0.5, 0, 0.5],
|
|
[0.5, 0.5, 0]],
|
|
supercell_matrix=np.diag([2, 2, 2]),
|
|
symprec=1e-5)
|
|
|
|
nac_params = {'born': borns_,
|
|
'factor': Hartree * Bohr,
|
|
'dielectric': epsilon_}
|
|
```
|