.. |umulaut_u| unicode:: U+00FC ANPHON: Input files ------------------- Format of input files ~~~~~~~~~~~~~~~~~~~~~ Each input file should consist of entry fields. Available entry fields are **&general**, **&cell**, **&scph**, **&qha**, **&relax**, **&kpoint**, **&strain**, and **&displace**. The format of the input file is the same as that of *alm* which can be found :ref:`here `. .. _label_inputvar_anphon: List of supported input variables ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. csv-table:: :widths: 25, 25, 25, 25 **&general** :ref:`BCONNECT `, :ref:`BORNINFO `, :ref:`BORNSYM `, :ref:`CLASSICAL ` :ref:`EMIN `, :ref:`EPSILON `, :ref:`FC2XML `, :ref:`FCSXML ` :ref:`ISMEAR `, :ref:`KD `, :ref:`MASS `, :ref:`MODE ` :ref:`NA_SIGMA `, :ref:`NKD `, :ref:`NONANALYTIC `, :ref:`PREFIX ` :ref:`PRINTSYM `, :ref:`RESTART `, :ref:`TMIN `, :ref:`TOLERANCE ` :ref:`TRISYM ` **&scph** :ref:`BUBBLE `, :ref:`IALGO `, :ref:`KMESH_INTERPOLATE `, :ref:`KMESH_SCPH ` :ref:`LOWER_TEMP `, :ref:`MAXITER `, :ref:`MIXALPHA `, :ref:`RELAX_STR ` :ref:`RESTART_SCPH `, :ref:`SELF_OFFDIAG `, :ref:`TOL_SCPH `, :ref:`WARMSTART ` **&qha** :ref:`KMESH_INTERPOLATE `, :ref:`KMESH_QHA `, :ref:`LOWER_TEMP `, :ref:`QHA_SCHEME ` :ref:`RELAX_STR ` **&relax** :ref:`ADD_HESS_DIAG `, :ref:`ALPHA_STDECENT `, :ref:`CELL_CONV_TOL `, :ref:`COOLING_U0_INDEX ` :ref:`COOLING_U0_THR `, :ref:`COORD_CONV_TOL `, :ref:`MAX_STR_ITER `, :ref:`MIXBETA_CELL ` :ref:`MIXBETA_COORD `, :ref:`RELAX_ALGO `, :ref:`RENORM_2TO1ST `, :ref:`RENORM_34TO1ST ` :ref:`RENORM_3TO2ND `, :ref:`SET_INIT_STR `, :ref:`STAT_PRESSURE `, :ref:`STRAIN_IFC_DIR ` **&analysis** :ref:`ANIME `, :ref:`ANIME_FRAMES `, :ref:`ANIME_CELLSIZE `, :ref:`GRUNEISEN ` :ref:`ISOFACT `, :ref:`ISOTOPE `, :ref:`KAPPA_COHERENT `, :ref:`KAPPA_SPEC ` :ref:`PDOS `, :ref:`PRINTEVEC `, :ref:`PRINTMSD `, :ref:`PRINTPR ` :ref:`PRINTVEL `, :ref:`PRINTXSF `, :ref:`SPS `, :ref:`TDOS ` :ref:`UCORR `, :ref:`ZMODE ` Description of input variables ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ "&general"-field ++++++++++++++++ .. _anphon_prefix: * **PREFIX**-tag : Job prefix to be used for names of output files :Default: None :Type: String ```` .. _anphon_mode: * **MODE**-tag = phonons | RTA ========= ============================================================== phonons | Calculate phonon dispersion relation, phonon DOS, | Gr\ |umulaut_u|\ neisen parameters etc. RTA | Calculate phonon lifetimes and lattice thermal conductivity | based on the Boltzmann transport equation (BTE) | with the relaxation time approximation (RTA). SCPH | Calculate temperature dependent phonon dispersion curves | by the self-consistent phonon method. ========= ============================================================== :Default: None :Type: String ```` .. _anphon_nkd: * **NKD**-tag : Number of atomic species :Default: None :Type: Integer ```` .. _anphon_kd: * **KD**-tag = Name[1], ... , Name[``NKD``] :Default: None :Type: Array of strings :Example: In the case of GaAs with ``NKD = 2``, it should be ``KD = Ga As``. ```` .. _anphon_mass: * MASS-tag = mass[1], ... , mass[``NKD``] :Default: Standard atomic weight of elements given by the ``KD``-tag :Type: Array of double :Example: In the case of Bi\ :sub:`2`\ Te\ :sub:`3` with ``NKD = 2``, ``MASS`` should be ``MASS = 208.98 127.60``. ```` .. _anphon_fcsxml: * **FCSXML**-tag : XML file containing force constants generated by the program *alm* :Default: None :Type: String ```` .. _anphon_fc2xml: * FC2XML-tag : XML file containing harmonic force constants for different size of supercell :Default: None :Type: String :Description: When ``FC2XML`` is given, the harmonic force constants in this file are used for calculating dynamical matrices. It is possible to use different size of supercell for harmonic and anharmonic terms, which are specified by ``FC2XML`` and ``FCSXML`` respectively. ```` .. _anphon_tolerance: * TOLERANCE-tag : Tolerance for finding symmetry operations :Default: 1.0e-6 :Type: Double ```` .. _anphon_printsym: * PRINTSYM-tag = 0 | 1 === ======================================================= 0 Symmetry operations won’t be saved in “SYMM_INFO_PRIM” 1 Symmetry operations will be saved in “SYMM_INFO_PRIM” === ======================================================= :Default: 0 :type: Integer ```` .. _anphon_nonanalytic: * NONANALYTIC-tag = 0 | 1 | 2 | 3 === =================================================================================== 0 | Non-analytic correction is not considered. 1 | Include the non-analytic correction by the damping method proposed by Parlinski. 2 | Include the non-analytic correction by the mixed-space approach 3 | Include the non-analytic correction by the Ewald method === =================================================================================== :Default: 0 :Type: Integer :Description: When ``NONANALYTIC > 0``, appropriate ``BORNINFO`` needs to be given. If ``NONANALYTIC = 1``, one may need to adjust the ``NA_SIGMA`` value to obtain reasonably smooth dispersion curves. ```` .. _anphon_na_sigma: * NA_SIGMA-tag : Damping factor for the non-analytic term :Default: 0.0 :Type: Double :Description: Used when ``NONANALYTIC = 1``. The definition of ``NA_SIGMA`` is described in the formalism section. ```` .. _anphon_borninfo: * BORNINFO-tag : File containing the macroscopic dielectric tensor and Born effective charges for the non-analytic correction :Default: None :Type: String :Description: The details of the file format can be found :ref:`here `. ```` .. _anphon_bornsym: * BORNSYM-tag = 0 | 1 === ================================================================= 0 Do not symmetrize Born effective charges 1 Symmetrize Born effective charges by using point group symmetry === ================================================================= :Default: 0 :Type: Integer ```` .. _anphon_tmin: * TMIN, TMAX, DT-tags : Temperature range and its stride in units of Kelvin :Default: ``TMIN = 0``, ``TMAX = 1000``, ``DT = 10`` :Type: Double ```` .. _anphon_emin: * EMIN, EMAX, DELTA_E-tags : Energy range and its stride in units of kayser (cm\ :sup:`-1`) :Default: ``EMIN`` and ``EMAX`` are set automatically from the eigenfrequencies as of ver. 1.5.0. The default value for ``DELTA_E`` is 10.0. :Type: Double ```` .. _anphon_ismear: * ISMEAR-tag = -1 | 0 | 1 === ======================================================= -1 Tetrahedron method 0 Lorentzian smearing with width of ``EPSILON`` 1 Gaussian smearing with width of ``EPSILON`` === ======================================================= :Default: -1 :Type: Integer :Description: ``ISMEAR`` specifies the method for Brillouin zone integration ```` .. _anphon_epsilon: * EPSILON-tag : Smearing width in units of Kayser (cm\ :sup:`-1`) :Default: 10.0 :Type: Double :Description: This variable is neglected when ``ISMEAR = -1`` ```` .. _anphon_bconnect: * BCONNECT-tag = 0 | 1 | 2 === =================================================================================== 0 | Phonon band is saved without change (sorted in order of energy) 1 | Phonon band is connected by using the similarity of eigenvectors. 2 | Same as ``BCONNECT=1``. In addition, information about the connectivity is | saved as ``PREFIX.connection``. === =================================================================================== :Default: 0 :Type: Integer :Description: The algorithm for connecting a band structure is described here_. .. _here : https://www.slideshare.net/TakeshiNishimatsu/two-efficient-algorithms-for-drawing-accurate-and-beautiful-phonon-dispersion ```` .. _anphon_classical: * CLASSICAL-tag = 0 | 1 === ======================================================= 0 Use quantum statistics (default) 1 Use classical statistics === ======================================================= :Default: 0 :Type: Integer :Description: When ``CLASSICAL = 1``, all thermodynamic functions including the occupation function, heat capacity, and mean square displacements are calculated using the classical formulae. This option may be useful when comparing the lattice dynamics and molecular dynamics results. .. list-table:: Comparison of quantum and classical values :header-rows: 1 * - Function - Quantum (``CLASSICAL = 0``) - Classical (``CLASSICAL = 1``) * - Occupation number - :math:`\displaystyle n_\mathrm{B}=\frac{1}{\exp(\beta\hbar\omega) - 1}` - :math:`\displaystyle n_\mathrm{C}=\frac{1}{\beta\hbar\omega}` * - Mode specific heat - :math:`\displaystyle c_{q} = k_{\mathrm{B}}\left[\frac{\beta\hbar\omega_q}{2}\mathrm{csch}\bigg({\frac{\beta\hbar\omega_q}{2}}\bigg)\right]^2` - :math:`\displaystyle c_{q} = k_{\mathrm{B}}` * - MSD of normal mode :math:`\braket{Q^{*}_qQ_q}` - :math:`\displaystyle \frac{\hbar (1 + n_{\mathrm{B}})}{2\omega_q}` - :math:`\displaystyle \frac{1}{\beta\omega_{q}^{2}}` ```` .. _anphon_trisym: * TRISYM-tag : Flag to use symmetry operations to reduce the number of triples of :math:`k` points for self-energy calculations === ======================================================= 0 Symmetry will not be used 1 Use symmetry to reduce triples of :math:`k` points === ======================================================= :Default: 1 :Type: Integer :Description: This variable is used only when ``MODE = RTA``. .. Note:: ``TRISYM = 1`` can reduce the computational cost, but phonon linewidth stored to the file ``PREFIX``.result needs to be averaged at points of degeneracy. For that purpose, a subsidiary program *analyze_phonons.py** should be used. ```` .. _anphon_restart: * RESTART-tag : Flag to restart the calculation when ``MODE = RTA`` === ======================================================= 0 Calculate from scratch 1 Restart from the existing file === ======================================================= :Default: 1 if there is a file named ``PREFIX``.result; 0 otherwise :Type: Integer ```` "&scph"-field (Read only when ``MODE = SCPH``) ++++++++++++++++++++++++++++++++++++++++++++++ .. _anphon_kmesh_interpolate: * KMESH_INTERPOLATE-tag = k1, k2, k3 :Default: None :Type: Array of integers :Description: In the iteration process of the SCPH equation, the interpolation is done using the :math:`k` mesh defined by ``KMESH_INTERPOLATE``. ```` .. _anphon_kmesh_scph: * KMESH_SCPH-tag = k1, k2, k3 :Default: None :Type: Array of integers :Description: This :math:`k` mesh is used for the inner loop of the SCPH equation. Each value of ``KMESH_SCPH`` must be equal to or a multiple of the number of ``KMESH_INTERPOLATE`` in the same direction. ```` .. _anphon_self_offdiag: * SELF_OFFDIAG-tag = 0 | 1 === ================================================================================ 0 Neglect the off-diagonal elements of the loop diagram in the SCPH calculation 1 Consider the off-diagonal elements of the loop diagram in the SCPH calculation === ================================================================================ :Default: 0 :Type: Integer :Description: ``SELF_OFFDIAG = 1`` is more accurate, but expensive. ```` .. _anphon_tol_scph: * TOL_SCPH-tag: Stopping criterion of the SCPH iteration :Default: 1.0e-10 :Type: Double :Description: The SCPH iteration stops when both :math:`[\frac{1}{N_{q}}\sum_{q} (\Omega_{q}^{(i)}-\Omega_{q}^{(i-1)})^{2}]^{1/2}` < ``TOL_SCPH`` and :math:`(\Omega_{q}^{(i)})^{2} \geq 0 \; (\forall q)` are satisfied. Here, :math:`\Omega_{q}^{(i)}` is the anharmonic phonon frequency in the :math:`i`\ th iteration and :math:`q` is the phonon modes at the irreducible momentum grid of ``KMESH_INTERPOLATE``. ```` .. _anphon_mixalpha: * MIXALPHA-tag: Mixing parameter used in the SCPH iteration :Default: 0.1 :Type: Double ```` .. _anphon_maxiter: * MAXITER-tag: Maximum number of the SCPH iteration :Default: 1000 :Type: Integer ```` .. _anphon_lower_temp: * LOWER_TEMP-tag = 0 | 1 === =============================================================================== 0 The SCPH iteration start from ``TMIN`` to ``TMAX``. (Raise the temperature) 1 The SCPH iteration start from ``TMAX`` to ``TMIN``. (Lower the temperature) === =============================================================================== :Default: 1 :Type: Integer ```` .. _anphon_warmstart: * WARMSTART-tag = 0 | 1 === =============================================================================== 0 SCPH iteration is initialized by harmonic frequencies and eigenvectors 1 SCPH iteration is initialized by the solution of the previous temperature === =============================================================================== :Default: 1 :Type: Integer :Description: ``WARMSTART = 1`` usually improves the convergence. ```` .. _anphon_ialgo: * IALGO-tag = 0 | 1 === =============================================================================== 0 MPI parallelization for the :math:`k` point 1 MPI parallelization for the phonon branch === =============================================================================== :Default: 0 :Type: Integer :Description: Use ``IALGO = 1`` when the primitive cell contains many atoms and the number of :math:`k` points is small. ```` .. _anphon_restart_scph: * RESTART_SCPH-tag = 0 | 1 === ============================================================== 0 Perform a SCPH calculation from scratch 1 Skip a SCPH iteration by loading a precalculated result === ============================================================== :Default: 1 if the file ``PREFIX.scph_dymat`` exists in the working directory; 0 otherwise :Type: Integer ```` .. _anphon_bubble: * BUBBLE-tag = 0 | 1 === ============================================================== 0 No bubble correction to the dynamical matrix 1 Calculate bubble correction on top of the SCPH dynamical matrix === ============================================================== :Default: 0 :Type: Integer ```` .. _anphon_relax_str: * RELAX_STR-tag = 0 | 1 | 2 | 3 === ============================================================== 0 Don't relax the crystal structure (not supported when ``mode = QHA``). 1 Relax atomic positions. 2 Relax both atomic positions and the shape of the unit cell. 3 Lowest-order perturbation theory (not supported when ``MODE = SCPH``). === ============================================================== :Default: 0 :Type: Integer ```` "&qha"-field (Read only when ``MODE = QHA``) ++++++++++++++++++++++++++++++++++++++++++++++ .. _anphon_qha_kmesh_interpolate: * KMESH_INTERPOLATE-tag = k1, k2, k3 :Default: None :Type: Array of integers :Description: In the structural optimization based on quasiharmonic approximation (QHA), the interpolation is done using the :math:`k` mesh defined by ``KMESH_INTERPOLATE``. ```` .. _anphon_qha_kmesh_qha: * KMESH_QHA-tag = k1, k2, k3 :Default: None :Type: Array of integers :Description: This :math:`k` mesh is used for the QHA-based structural optimization. Each value of ``KMESH_QHA`` must be equal to or a multiple of the number of ``KMESH_INTERPOLATE`` in the same direction. ```` .. _anphon_qha_relax_str: * RELAX_STR-tag = 0 | 1 | 2 | 3 === ============================================================== 0 Don't relax the crystal structure (not supported when ``mode = QHA``). 1 Relax atomic positions. 2 Relax both atomic positions and the shape of the unit cell. 3 Lowest-order perturbation theory (not supported when ``mode = SCPH``). === ============================================================== :Default: 0 :Type: Integer ```` .. _anphon_qha_lower_temp: * LOWER_TEMP-tag = 0 | 1 === =============================================================================== 0 The structural optimization start from ``TMIN`` to ``TMAX``. (Raise the temperature) 1 The structural optimization start from ``TMAX`` to ``TMIN``. (Lower the temperature) === =============================================================================== :Default: 1 :Type: Integer ```` .. _anphon_qha_scheme: * QHA_SCHEME-tag = 0 | 1 | 2 === ============================================================== 0 Full optimization within QHA. 1 zero-static internal stress approximation (ZSISA). 2 volumetric ZSISA (v-ZSISA). === ============================================================== :Default: 0 :Type: Integer :Description: This option is used only when ``mode = QHA`` and ``RELAX_STR = 2``. ```` "&relax"-field (Read only when ``RELAX_STR != 0``) ++++++++++++++++ .. _anphon_relax_algo: * RELAX_ALGO-tag = 1 | 2 === ============================================================== 1 Steepest decent (not recommended) 2 Newton-like method === ============================================================== :Default: 2 :Type: Integer :Description: Algorithm to update the crystal structure in structural optimization. This option is used only when ``RELAX_STR = 1, 2``. ``RELAX_ALGO = 1`` works properly only when the unit cell is fixed (``RELAX_STR = 1``). ```` .. _anphon_alpha_stdecent: * ALPHA_STDECENT-tag: Coefficient of steepest decent in structural optimization :Default: 1.0e4 :Type: Double :Description: :math:`\alpha` coefficient in structural optimization with steepest-decent algorithm. The unit is [:math:`m_e a_B^2/(2\text{Ry})`]. This option is used only when ``RELAX_ALGO = 1``. ```` .. _anphon_max_str_iter: * MAX_STR_ITER-tag: Maximum number of structure updates. :Default: 100 :Type: Integer :Description: This option is used only when ``RELAX_STR = 1, 2``. ```` .. _anphon_add_hess_diag: * ADD_HESS_DIAG-tag: Correction to the estimated Hessian of free energy in units of kayser (cm\ :sup:`-1`) :Default: 100.0 :Type: Double :Description: The squared ``ADD_HESS_DIAG`` is added to the diagonal components of estimated Hessians, which is used to update crystal structures in structural optimization. ``ADD_HESS_DIAG`` makes the calculation more robust in the presence of soft modes near the structural phase transition, but setting large values will make the convergence slower. This option is used only when ``RELAX_ALGO = 2``. ```` .. _anphon_coord_conv_tol: * COORD_CONV_TOL-tag: Threshold of convergence for atomic positions in structural optimization. :Default: 1.0e-5 :Type: Double :Description: The value is interpreted in units of Bohr. This option is used only when ``RELAX_STR = 1, 2``. ```` .. _anphon_mixbeta_coord: * MIXBETA_COORD-tag: Mixing coefficient for atomic positions in structure updates. :Default: 0.5 :Type: Double :Description: This option is used only when ``RELAX_STR = 1, 2``. ```` .. _anphon_cell_conv_tol: * CELL_CONV_TOL-tag: Threshold of convergence for displacement gradient tensor :math:`u_{\mu \nu}` in structural optimization. :Default: 1.0e-5 :Type: Double :Description: This option is used only when ``RELAX_STR = 2``. ```` .. _anphon_mixbeta_cell: * MIXBETA_CELL-tag: Mixing coefficient for displacement gradient tensor :math:`u_{\mu \nu}` in structure updates. :Default: 0.5 :Type: Double :Description: This option is used only when ``RELAX_STR = 2``. ```` .. _anphon_set_init_str: * SET_INIT_STR-tag = 1 | 2 | 3 === ============================================================== 1 Set initial structure from the input file at each temperature. 2 Start from the crystal structure of the previous temperature. 3 Start from the crystal structure of the previous temperature in low-symmetry phase. === ============================================================== :Default: 1 :Type: Integer :Description: This option specifies how to set the initial structure of structural optimization at different temperatures. This option is used when ``RELAX_STR = 1, 2``. In all options, the initial structure at the initial temperature is set from the input file. The initial structure of the input file is read from the ``&strain`` and ``&displace`` field. When ``SET_INIT_STR = 3``, the initial displacement from the input file is used if the crystal structure converges to the high-symmetry phase in the previous temperature. The criteria to distinguish low-symmetry and high-symmetry phases is explained in :ref:`COOLING_U0_THR `. ```` .. _anphon_cooling_u0_index: * COOLING_U0_INDEX-tag = 0 | 1 | ... | 3N-1 (N : the number of atoms in the unit cell) :Default: 0 :Type: Integer :Description: Specify as :math:`3\times\alpha + \mu`. Here, :math:`\alpha` denotes the atom index in the primitive cell and :math:`\mu` is the xyz index, where both indices are zero-indexed. See the description of :ref:`COOLING_U0_THR ` for details. This option is used only when ``SET_INIT_STR = 3``. ```` .. _anphon_cooling_u0_thr: * COOLING_U0_THR-tag: Threshold to judge high-symmetry phase in structural optimization [Bohr]. :Default: 0.001 :Type: Double :Description: The crystal structure is judged to be back to the high-symmetry phase if :math:`u^{(0)}` [``COOLING_U0_INDEX``] < ``COOLING_U0_THR``. This option is useful in cooling calculations because small displacements to the high-symmetry structure is required to induce spontaneous symmetry breaking. This option is used only when ``SET_INIT_STR = 3``. ```` .. _anphon_stat_pressure: * STAT_PRESSURE-tag: Hydrostatic pressure in GPa. :Default: 0.0 :Type: Double ```` .. _anphon_renorm_2to1st: * RENORM_2TO1ST-tag = 0 | 1 | 2 === ============================================================== 0 Set zero. 1 Real-space IFC renormalization. (not recommended) 2 Finite difference method with respect to strain. === ============================================================== :Default: 2 :Type: Integer :Description: This option specifies the method to calculate first-order derivatives of first-order IFCs with respect to strain :math:`\frac{\partial \Phi_{\mu}(0\alpha)}{\partial u_{\mu_1 \nu_1} }`. This option is used only when ``RELAX_STR = 2, 3``. Note that ``RENORM_2TO1ST = 1`` requires rotational invariance on IFCs, which is not checked in the program ANPHON. ``RENORM_2TO1ST = 0`` can be used for high-symmetry materials in which strain-force coupling is zero, which a user need to confirm themselves. ```` .. _anphon_renorm_34to1st: * RENORM_34TO1ST-tag = 0 | 1 === ============================================================== 0 Set zero. 1 Real-space IFC renormalization. === ============================================================== :Default: 0 :Type: Integer :Description: This option specifies the method to calculate second and higher-order derivatives of first-order IFCs with respect to strain. :math:`\frac{\partial^2 \Phi_{\mu}(0\alpha)}{\partial u_{\mu_1 \nu_1} \partial u_{\mu_2 \nu_2}}`, :math:`\frac{\partial^3 \Phi_{\mu}(0\alpha)}{\partial u_{\mu_1 \nu_1} \partial u_{\mu_2 \nu_2} \partial u_{\mu_3 \nu_3}}` This option is used only when ``RELAX_STR = 2, 3``. Note that ``RENORM_34TO1ST = 1`` requires rotational invariance on IFCs, which a user need to confirm themselves. ```` .. _anphon_renorm_3to2nd: * RENORM_3TO2ND-tag = 1 | 2 | 3 === ============================================================== 1 Real-space IFC renormalization. 2 Finite difference method (Read input from all six strain patterns). 3 Finite difference method (Read input from specified strain patterns). === ============================================================== :Default: 2 :Type: Integer :Description: This option specifies the method to calculate first-order derivatives of harmonic IFCs with respect to strain. :math:`\frac{\partial \Phi_{\mu_1 \mu_2}(0\alpha_1, R \alpha_2)}{\partial u_{\mu \nu}}` This option is used only when ``RELAX_STR = 2, 3``. To use ``RENORM_3TO2ND = 3``, the entries of the rotation matrices of ALL the space-group operations must be either 0 or :math:`\pm` 1 in Cartesian representation. ```` .. _anphon_strain_ifc_dir: * STRAIN_IFC_DIR-tag: Directory name of the inputs of strain-IFC couplings. :Default: None :Type: String :Description: When ``RENORM_2TO1ST = 2 `` or ``RENORM_3TO2ND = 3``, the input files of the strain-IFC couplings must be given properly in this directory. ```` "&cell"-field +++++++++++++ Please specify the cell parameters of the *primitive cell* as:: &cell a a11 a12 a13 a21 a22 a23 a31 a32 a33 / The cell parameters are then given by :math:`\vec{a}_{1} = a \times (a_{11}, a_{12}, a_{13})`, :math:`\vec{a}_{2} = a \times (a_{21}, a_{22}, a_{23})`, and :math:`\vec{a}_{3} = a \times (a_{31}, a_{32}, a_{33})`. .. Note:: The lattice constant :math:`a` must be consistent with the value used for the program *alm*. For example, if one used :math:`a = 20.4 a_{0}` for a 2x2x2 supercell of Si, one should use :math:`a = 10.2 a_{0}` here for the primitive cell. ```` "&kpoint"-field +++++++++++++++ This entry field is used to specify the list of :math:`k` points to be calculated. The first entry **KPMODE** specifies the types of calculation which is followed by detailed entries. * **KPMODE = 0** : Calculate phonon frequencies at given :math:`k` points For example, if one wants to calculate phonon frequencies at Gamma (0, 0, 0) and X (0, 1/2, 1/2) of an FCC crystal, the ``&kpoint`` entry should be written as :: &kpoint 0 0.000 0.000 0.000 0.000 0.500 0.500 / * **KPMODE = 1** : Band dispersion calculation For example, if one wants to calculate phonon dispersion relations along G\-K\-X\-G\-L of a FCC crystal, the ``&kpoint`` entry should be written as follows:: &kpoint 1 G 0.000 0.000 0.000 K 0.375 0.375 0.750 51 K 0.375 0.375 0.750 X 0.500 0.500 1.000 51 X 0.000 0.500 0.500 G 0.000 0.000 0.000 51 G 0.000 0.000 0.000 L 0.500 0.500 0.500 51 / The 1st and 5th columns specify the character of Brillouin zone edges, which are followed by fractional coordinates of each point. The last column indicates the number of sampling points. * **KPMODE = 2** : Uniform :math:`k` grid for phonon DOS and thermal conductivity In order to perform a calculation with 20x20x20 :math:`k` grid, the entry should be :: &kpoint 2 20 20 20 / ```` "&strain"-field (Read only when ``RELAX_STR = 2``) +++++++++++++++++ Please specify the initial displacement gradient tensor :math:`u_{\mu \nu}` for structural optimization as :: &cell u_xx u_xy u_xz u_yx u_yy u_yz u_zx u_zy u_zz / Note that a user need to give a symmetric matrix. "&displace"-field (Read only when ``RELAX_STR = 1, 2``) +++++++++++++++++ Please specify the initial atomic displacements :math:`u^{(0)}_{\alpha \mu}` [Bohr]. * **DISPMODE = 0** : Fractional coordinate representation The ``&displace`` entry should be written as follows. The first four lines after DISPMODE (= 0) specifies the unit cell, whose format is the same as the ``&cell`` field. Note that the unit cell in the ``&displace`` field is used only for transforming the input to the real space representation. Thus, the unit cell here does not need to be commensurate with the primitive cell or some supercells. u_ij is the j-th component of the displacement of i-th atom in the primitive cell in fractional coordinate representation. :: &displace 0 a a11 a12 a13 a21 a22 a23 a31 a32 a33 u_01, u_02, u_03 ... / * **DISPMODE = 1** : Cartesian coordinate representation Each line after DISPMODE (= 1) specifies the initial atomic displacement in Cartesian representation. u_ij is the j component of the displacement of i-th atom in the primitive cell. :: &displace 1 u_0x, u_0y, u_0z ... / "&analysis"-field +++++++++++++++++ .. _anphon_gruneisen: * GRUNEISEN-tag = 0 | 1 === =================================================================== 0 Gr\ |umulaut_u|\ neisen parameters will not be calculated 1 Gr\ |umulaut_u|\ neisen parameters will be stored === =================================================================== :Default: 0 :Type: Integer :Description: When ``MODE = phonons`` and ``GRUNEISEN = 1``, Gr\ |umulaut_u|\ neisen parameters will be stored in ``PREFIX``.gru (*KPMODE* = 1) or ``PREFIX``.gru_all (*KPMODE* = 2). .. Note:: To compute Gr\ |umulaut_u|\ neisen parameters, cubic force constants must be contained in the ``FCSXML`` file. ```` .. _anphon_printevec: * PRINTEVEC-tag = 0 | 1 === =================================================================== 0 Do not print phonon eigenvectors 1 Print phonon eigenvectors in the ``PREFIX``.evec file === =================================================================== :Default: 0 :Type: Integer ```` .. _anphon_printxsf: * PRINTXSF-tag = 0 | 1 === =================================================================== 0 Do not save an AXSF file 1 Create an AXSF file ``PREFIX``.axsf === =================================================================== :Default: 0 :Type: Integer :Description: This is to visualize the direction of vibrational modes at gamma (0, 0, 0) by XCrySDen. This option is valid only when ``MODE = phonons``. ```` .. _anphon_printvel: * PRINTVEL-tag = 0 | 1 === =================================================================== 0 Do not print group velocity 1 Store phonon velocities to a file === =================================================================== :Default: 0 :Type: Integer :Description: When ``MODE = phonons`` and ``PRINTVEL = 1``, group velocities of phonons will be stored in ``PREFIX``.phvel (*KPMODE* = 1) or ``PREFIX``.phvel_all (*KPMODE* = 2). ```` .. _anphon_printmsd: * PRINTMSD-tag = 0 | 1 === =================================================================== 0 Do not print mean-square-displacement (MSD) of atoms 1 Save MSD of atoms to the file ``PREFIX``.mds === =================================================================== :Default: 0 :Type: Integer :Description: This flag is available only when ``MODE = phonons`` and *KPMODE* = 2. ```` .. _anphon_pdos: * PDOS-tag = 0 | 1 === =================================================================== 0 Only the total DOS will be printed in ``PREFIX``.dos 1 Atom-projected phonon DOS will be stored in ``PREFIX``.dos === =================================================================== :Default: 0 :Type: Integer :Description: This flag is available only when ``MODE = phonons`` and *KPMODE* = 2. ```` .. _anphon_tdos: * TDOS-tag = 0 | 1 === =================================================================== 0 Do not compute two-phonon DOS 1 Two-phonon DOSs will be stored in ``PREFIX``.tdos === =================================================================== :Default: 0 :Type: Integer :Description: This flag is available only when ``MODE = phonons`` and *KPMODE* = 2. .. Note:: Calculation of two-phonon DOS is computationally expensive. ```` .. _anphon_sps: * SPS-tag = 0 | 1 | 2 === ==================================================================================== 0 Do not compute scattering phase space 1 | Total and mode-decomposed scattering phase space involving | the three-phonon processes will be stored in ``PREFIX``.sps 2 Three-phonon scattering phase space with the Bose factor will be stored in ``PREFIX``.sps_Bose === ==================================================================================== :Default: 0 :Type: Integer :Description: This flag is available only when ``MODE = phonons`` and *KPMODE* = 2. ```` .. _anphon_printpr: * PRINTPR-tag = 0 | 1 === ==================================================================================== 0 Do not compute the (atomic) participation ratio 1 | Compute participation ratio and atomic participation ratio, which will be | stored in ``PREFIX``.pr and ``PREFIX``.apr respectively. === ==================================================================================== :Default: 0 :Type: Integer :Description: This flag is available when ``MODE = phonons``. ```` .. _anphon_kappa_coherent: * KAPPA_COHERENT-tag = 0 | 1 | 2 === ==================================================================================== 0 Do not compute the coherent component of thermal conductivity 1 Compute the coherent component of thermal conductivity and save it in ``PREFIX``.kl_coherent. 2 | In addition to above (``KAPPA_COHERENT = 1``), all elements of the coherent term | are saved in ``PREFIX``.kc_elem. === ==================================================================================== :Default: 0 :Type: Integer :Description: This flag is available when ``MODE = RTA``. For the theoretical details, please see :ref:`this page `. .. caution:: Still experimental. Please check the validity of results carefully. ```` .. _anphon_kappa_spec: * KAPPA_SPEC-tag = 0 | 1 === ==================================================================================== 0 Do not compute the thermal conductivity spectra 1 Compute the thermal conductivity spectra, which will be stored in ``PREFIX``.kappa_spec . === ==================================================================================== :Default: 0 :Type: Integer :Description: This flag is available when ``MODE = RTA``. ```` .. _anphon_isotope: * ISOTOPE-tag = 0 | 1 === ========================================================================= 0 Do not consider phonon-isotope scatterings 1 Consider phonon-isotope scatterings 2 | Consider phonon-isotope scatterings as in ``ISOTOPE = 1`` and | the calculated selfenergy is stored in ``PREFIX``.gamma_isotope === ========================================================================= :Default: 0 :Type: Integer :Description: When ``MODE = RTA`` and ``ISOTOPE = 1 or 2``, phonon scatterings due to isotopes will be considered perturbatively. ``ISOFACT`` should be properly given. ```` .. _anphon_isofact: * ISOFACT-tag = isofact[1], ... , isofact[``NKD``] :Default: Automatically calculated from the ``KD`` tag :Type: Array of doubles :Description: Isotope factor is a dimensionless value defined by :math:`\sum_{i} f_{i} (1 - m_{i}/\bar{m})^{2}`. Here, :math:`f_{i}` is the fraction of the :math:`i`\ th isotope of an element having mass :math:`m_{i}`, and :math:`\bar{m}=\sum_{i}f_{i}m_{i}` is the average mass, respectively. This quantity is equivalent to :math:`g_{2}` appearing in the original paper by S. Tamura [Phys. Rev. B, 27, 858.]. ```` .. _anphon_ucorr: * UCORR-tag = 0 | 1 === ========================================================================= 0 Do nothing 1 | Compute the displacement-displacement correlation function. | The result is stored in ``PREFIX``.ucorr === ========================================================================= :Default: 0 :Type: Integer :Description: The displacement-displacement correlation function involves two atoms. The first atom is located in the primitive cell at the center (shift1=[0,0,0]) and the second atom is located in the :math:`\ell'`\ th cell. The translation vector to the :math:`\ell'`\ th cell can be specified by the ``SHIFT_UCORR`` tag. This tag is effective only when ``MODE = phonons`` and *KPMODE* = 2 ```` .. _anphon_shift_ucorr: * SHIFT_UCORR-tag = l1, l2, l3 :Default: [0, 0, 0] :Type: Array of integers :Description: This tag specifies the translation vector used for computing the displacement-displacement (uu) correlation function. For example, if one wants to compute the uu correlation function between an atom 1 in the cell at the center and atom 2 in the neighboring cell at :math:`\boldsymbol{r}(\ell')=(1,0,0)`, ``SHIFT_UCORR`` should be set as ``SHIFT_UCORR = 1 0 0``. ```` .. _anphon_zmode: * ZMODE-tag = 0 | 1 === ========================================================================= 0 Do nothing 1 | Compute the mode effective charges of the zone-center phonons. | The result is stored in ``PREFIX``.zmode === ========================================================================= :Default: 0 :Type: Integer :Description: When ``MODE = phonons`` and ``ZMODE = 1``, the mode effective charges are computed for the phonon modes at the Gamma point and saved in ``PREFIX``.zmode. The unit of the mode effective charge is :math:`e \; \text{amu}^{-1/2}`. ```` .. .. _anphon_fe_bubble: .. * FE_BUBBLE-tag = 0 | 1 .. === ==================================================================================== .. 0 Do not compute the vibrational free-energy associated with the bubble diagram .. 1 | Compute the vibrational free-energy associated with the bubble diagram and .. | save it in ``PREFIX``.thermo (when ``MODE = phonons``) or ``PREFIX``.scph_thermo (when ``MODE = SCPH``). .. === ==================================================================================== .. :Default: 0 .. :Type: Integer .. :Description: This tag is used when *KPMODE* = 2. .. ```` .. _anphon_anime: * ANIME-tag = k1, k2, k3 :Default: None :Type: Array of doubles :Description: This tag is to animate vibrational mode. k1, k2, and k3 specify the momentum of phonon modes to animate, which should be given in units of the reciprocal lattice vector. For example, ``ANIME = 0.0 0.0 0.5`` will animate phonon modes at (0, 0, 1/2). When ``ANIME`` is given, ``ANIME_CELLSIZE`` is also necessary. You can choose the format of animation files, either AXSF or XYZ, by ``ANIME_FORMAT`` tag. ```` .. _anphon_anime_frames: * ANIME_FRAMES-tag: The number of frames saved in animation files :Default: 20 :Type: Integer ```` .. _anphon_anime_cellsize: * ANIME_CELLSIZE-tag = L1, L2, L3 :Default: None :Type: Array of integers :Description: This tag specifies the cell size for animation. L1, L2, and L3 should be large enough to be commensurate with the reciprocal point given by the ``ANIME`` tag. ```` .. _anphon_anime_format: * ANIME_FORMAT = xsf | xyz :Default: xyz :Type: String :Description: When ``ANIME_FORMAT = xsf``, ``PREFIX``.anime???.axsf files are created for XcrySDen. When ``ANIME_FORMAT = xyz``, ``PREFIX``.anime???.xyz files are created for VMD (and any other supporting software such as Jmol). ```` .. _label_format_BORNINFO: Format of BORNINFO ~~~~~~~~~~~~~~~~~~ When one wants to consider the LO-TO splitting near the :math:`\Gamma` point, it is necessary to set ``NONANALYTIC > 0`` and provide ``BORNINFO`` file containing dielectric tensor :math:`\epsilon^{\infty}` and Born effective charge :math:`Z^{*}`. In ``BORNINFO`` file, the dielectric tensor should be written in first 3 lines which are followed by Born effective charge tensors for each atom as the following. .. math:: :nowrap: \begin{eqnarray*} \epsilon_{xx}^{\infty} & \epsilon_{xy}^{\infty} & \epsilon_{xz}^{\infty} \\ \epsilon_{yx}^{\infty} & \epsilon_{yy}^{\infty} & \epsilon_{yz}^{\infty} \\ \epsilon_{zx}^{\infty} & \epsilon_{zy}^{\infty} & \epsilon_{zz}^{\infty} \\ Z_{1,xx}^{*} & Z_{1,xy}^{*} & Z_{1,xz}^{*} \\ Z_{1,yx}^{*} & Z_{1,yy}^{*} & Z_{1,zz}^{*} \\ Z_{1,zx}^{*} & Z_{1,zy}^{*} & Z_{1,zz}^{*} \\ & \vdots & \\ Z_{N_p,xx}^{*} & Z_{N_p,xy}^{*} & Z_{N_p,xz}^{*} \\ Z_{N_p,yx}^{*} & Z_{N_p,yy}^{*} & Z_{N_p,zz}^{*} \\ Z_{N_p,zx}^{*} & Z_{N_p,zy}^{*} & Z_{N_p,zz}^{*} \\ \end{eqnarray*} Here, :math:`N_p` is the number of atoms contained in the *primitive cell*. .. Attention:: Please pay attention to the order of Born effective charges.