2012-12-11 16:19:23 +08:00
|
|
|
#!/usr/bin/env python
|
|
|
|
|
|
|
|
# Copyright (C) 2011 Atsushi Togo
|
|
|
|
# All rights reserved.
|
|
|
|
#
|
|
|
|
# This file is part of phonopy.
|
|
|
|
#
|
|
|
|
# Redistribution and use in source and binary forms, with or without
|
|
|
|
# modification, are permitted provided that the following conditions
|
|
|
|
# are met:
|
|
|
|
#
|
|
|
|
# * Redistributions of source code must retain the above copyright
|
|
|
|
# notice, this list of conditions and the following disclaimer.
|
|
|
|
#
|
|
|
|
# * Redistributions in binary form must reproduce the above copyright
|
|
|
|
# notice, this list of conditions and the following disclaimer in
|
|
|
|
# the documentation and/or other materials provided with the
|
|
|
|
# distribution.
|
|
|
|
#
|
|
|
|
# * Neither the name of the phonopy project nor the names of its
|
|
|
|
# contributors may be used to endorse or promote products derived
|
|
|
|
# from this software without specific prior written permission.
|
|
|
|
#
|
|
|
|
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
|
|
|
|
# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
|
|
|
|
# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
|
|
|
|
# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
|
|
|
|
# COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
|
|
|
# INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
|
|
|
|
# BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
|
|
|
|
# LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
|
|
|
|
# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
|
|
|
|
# LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
|
|
|
# ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
|
|
|
|
# POSSIBILITY OF SUCH DAMAGE.
|
|
|
|
|
|
|
|
import os
|
|
|
|
import sys
|
|
|
|
import numpy as np
|
|
|
|
from optparse import OptionParser
|
|
|
|
|
|
|
|
from phonopy.interface.vasp import read_vasp
|
|
|
|
from phonopy.structure.cells import get_supercell, Primitive, print_cell
|
|
|
|
from phonopy.structure.symmetry import Symmetry
|
|
|
|
from phonopy.harmonic.forces import Forces
|
|
|
|
from phonopy.harmonic.force_constants import get_force_constants, set_permutation_symmetry
|
|
|
|
from phonopy.harmonic.dynamical_matrix import DynamicalMatrix
|
|
|
|
from phonopy.hphonopy.file_IO import parse_FORCE_SETS
|
|
|
|
from phonopy.units import THzToCm, VaspToTHz
|
|
|
|
from phonopy.hphonopy.file_IO import parse_BORN
|
|
|
|
from anharmonic.fc3 import get_fc3, symmetrize_fc3
|
|
|
|
from anharmonic.file_IO import write_fc2_dat, write_fc3_dat,\
|
|
|
|
parse_disp_yaml, parse_FORCES_THIRD, write_FORCES_THIRD,\
|
|
|
|
parse_DELTA_FORCES, write_supercells_with_displacements,\
|
|
|
|
parse_QPOINTS3, parse_fc3, parse_fc2,\
|
|
|
|
write_DELTA_FC2_SETS, parse_DELTA_FC2_SETS,\
|
|
|
|
write_FC2_FOURTH_SETS, parse_FC2_FOURTH_SETS
|
|
|
|
from anharmonic.displacement_fc3 import get_third_order_displacements
|
|
|
|
from anharmonic.settings import Phono3pyConfParser
|
2013-02-18 14:05:15 +08:00
|
|
|
from anharmonic import Phono3py, JointDOS, get_gruneisen_parameters, \
|
|
|
|
get_ir_grid_points
|
2013-02-25 13:42:34 +08:00
|
|
|
from anharmonic.triplets import get_grid_address, reduce_grid_points
|
2012-12-11 16:19:23 +08:00
|
|
|
|
|
|
|
# AA is created at http://www.network-science.de/ascii/.
|
|
|
|
def print_phono3py():
|
|
|
|
print """ _ _____
|
|
|
|
_ __ | |__ ___ _ __ ___|___ / _ __ _ _
|
|
|
|
| '_ \| '_ \ / _ \| '_ \ / _ \ |_ \| '_ \| | | |
|
|
|
|
| |_) | | | | (_) | | | | (_) |__) | |_) | |_| |
|
|
|
|
| .__/|_| |_|\___/|_| |_|\___/____/| .__/ \__, |
|
|
|
|
|_| |_| |___/
|
|
|
|
"""
|
|
|
|
|
|
|
|
def print_end():
|
|
|
|
print """ _
|
|
|
|
___ _ __ __| |
|
|
|
|
/ _ \ '_ \ / _` |
|
|
|
|
| __/ | | | (_| |
|
|
|
|
\___|_| |_|\__,_|
|
|
|
|
"""
|
|
|
|
|
|
|
|
def print_error(message):
|
|
|
|
print message
|
|
|
|
|
|
|
|
# Parse options
|
|
|
|
parser = OptionParser()
|
|
|
|
parser.set_defaults(amplitude=None,
|
|
|
|
band_indices=None,
|
|
|
|
cell_poscar=None,
|
|
|
|
delta_fc2=False,
|
|
|
|
factor=None,
|
|
|
|
delta_fc2_sets_mode=False,
|
|
|
|
fc2_fourth_sets_mode=False,
|
|
|
|
freq_scale=None,
|
|
|
|
gamma_option=0,
|
|
|
|
grid_points=None,
|
|
|
|
is_cm=False,
|
|
|
|
is_decay_channel=False,
|
|
|
|
is_nodiag=False,
|
|
|
|
is_displacement=False,
|
|
|
|
is_nosym=False,
|
|
|
|
is_gruneisen=False,
|
|
|
|
is_joint_dos=False,
|
2013-02-08 17:34:46 +08:00
|
|
|
is_linewidth=False,
|
2013-02-18 14:05:15 +08:00
|
|
|
is_bterta=False,
|
2012-12-11 16:19:23 +08:00
|
|
|
is_nac=False,
|
|
|
|
is_plusminus_displacements=False,
|
|
|
|
is_read_triplets=False,
|
|
|
|
is_translational_symmetry=False,
|
|
|
|
is_symmetrize_fc2=False,
|
|
|
|
is_symmetrize_fc3_r=False,
|
|
|
|
is_symmetrize_fc3_q=False,
|
|
|
|
is_Peierls=False,
|
2013-02-18 14:05:15 +08:00
|
|
|
log_level=0,
|
2012-12-11 16:19:23 +08:00
|
|
|
mesh_numbers=None,
|
2013-02-25 13:42:34 +08:00
|
|
|
mesh_divisors=None,
|
2013-02-20 11:28:32 +08:00
|
|
|
multiple_sigmas=None,
|
|
|
|
no_kappa_stars=False,
|
2012-12-11 16:19:23 +08:00
|
|
|
q_direction=None,
|
|
|
|
primitive_axis=None,
|
|
|
|
read_fc2=False,
|
|
|
|
read_fc2_extra=False,
|
|
|
|
read_fc3=False,
|
2013-03-12 02:27:13 +08:00
|
|
|
read_gammas=False,
|
2012-12-11 16:19:23 +08:00
|
|
|
r2q_TI_index=None,
|
2013-02-20 17:21:35 +08:00
|
|
|
freq_step=None,
|
2012-12-11 16:19:23 +08:00
|
|
|
output_filename=None,
|
2013-02-05 17:00:25 +08:00
|
|
|
qpoints=None,
|
2012-12-11 16:19:23 +08:00
|
|
|
sigma=None,
|
|
|
|
supercell_dimension=None,
|
|
|
|
supercell_dimension_extra=None,
|
|
|
|
symprec=1e-5,
|
|
|
|
tmax=None,
|
|
|
|
tmin=None,
|
|
|
|
tstep=None,
|
2013-02-20 17:21:35 +08:00
|
|
|
temperatures=None,
|
2013-02-18 14:05:15 +08:00
|
|
|
verbose=True,
|
2013-03-02 17:52:22 +08:00
|
|
|
write_grid_points=False,
|
|
|
|
zheev_test=None)
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--amplitude", dest="amplitude", type="float",
|
|
|
|
help="Distance of displacements")
|
|
|
|
parser.add_option("--bi", "--band_indices", dest="band_indices",
|
|
|
|
type="string",
|
|
|
|
help="Band indices where life time is calculated")
|
2013-02-18 14:05:15 +08:00
|
|
|
parser.add_option("--br", "--bterta", dest="is_bterta",
|
|
|
|
action="store_true",
|
|
|
|
help="Calculate thermal conductivity in BTE-RTA")
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("-c", "--cell", dest="cell_poscar",
|
|
|
|
action="store", type="string",
|
|
|
|
help="Read unit cell", metavar="FILE")
|
|
|
|
parser.add_option("--cm", dest="is_cm",
|
|
|
|
action="store_true", help="Convert THz to cm^-1")
|
|
|
|
parser.add_option("-d", "--disp", dest="is_displacement",
|
|
|
|
action="store_true",
|
|
|
|
help="As first stage, get least displacements")
|
|
|
|
parser.add_option("--decay", dest="is_decay_channel",
|
|
|
|
action="store_true", help="Calculate decay channels")
|
|
|
|
parser.add_option("-o", dest="output_filename",
|
|
|
|
type="string",
|
|
|
|
help="Filename of output of damping function")
|
|
|
|
parser.add_option("--nodiag", dest="is_nodiag",
|
|
|
|
action="store_true",
|
|
|
|
help="Set displacements parallel to axes")
|
|
|
|
parser.add_option("--dim",
|
|
|
|
dest="supercell_dimension",
|
|
|
|
type="string",
|
|
|
|
help="Supercell dimension")
|
|
|
|
parser.add_option("--dim2",
|
|
|
|
dest="supercell_dimension_extra",
|
|
|
|
type="string",
|
|
|
|
help="Supercell dimension for extra fc2")
|
|
|
|
parser.add_option("--cf3", "--create_f3",
|
|
|
|
dest="forces_third_mode",
|
|
|
|
action="store_true",
|
|
|
|
help="Create FORCES_THIRD")
|
|
|
|
parser.add_option("--cdfc2", "--create_delta_fc2",
|
|
|
|
dest="delta_fc2_sets_mode",
|
|
|
|
action="store_true",
|
|
|
|
help="Create DELTA_FC2_SETS")
|
|
|
|
parser.add_option("--create_fc2_fourth",
|
|
|
|
dest="fc2_fourth_sets_mode",
|
|
|
|
action="store_true",
|
|
|
|
help="Create FC2_FOURTH_SETS")
|
|
|
|
parser.add_option("--factor", dest="factor", type="float",
|
|
|
|
help="Conversion factor to favorite frequency unit")
|
|
|
|
parser.add_option("--fc2",
|
|
|
|
dest="read_fc2",
|
|
|
|
action="store_true",
|
|
|
|
help="Read second order force constants")
|
|
|
|
parser.add_option("--fc2_extra",
|
|
|
|
dest="read_fc2_extra",
|
|
|
|
action="store_true",
|
|
|
|
help="Read extra second order force constants")
|
|
|
|
parser.add_option("--fc3",
|
|
|
|
dest="read_fc3",
|
|
|
|
action="store_true",
|
|
|
|
help="Read third order force constants")
|
|
|
|
parser.add_option("--delta_fc2",
|
|
|
|
dest="read_delta_fc2",
|
|
|
|
action="store_true",
|
|
|
|
help="Read DELTA_FC2_SETS")
|
|
|
|
parser.add_option("--fc2_fourth",
|
|
|
|
dest="read_fc2_fourth_sets",
|
|
|
|
action="store_true",
|
|
|
|
help="Read FC2_FOURTH_SETS")
|
|
|
|
parser.add_option("--freq_scale", dest="freq_scale", type="float",
|
2013-02-06 23:29:22 +08:00
|
|
|
help="Scale factor is multiplied to frequencies only, i.e., changes frequencies but assumed not to change the physical unit")
|
2013-02-18 14:05:15 +08:00
|
|
|
parser.add_option("--gamma_option", dest="gamma_option", type="int",
|
|
|
|
help="Option for the test of imaginary part of self energy")
|
|
|
|
parser.add_option("--gp", "--grid_points",
|
2012-12-11 16:19:23 +08:00
|
|
|
dest="grid_points",
|
|
|
|
type="string",
|
|
|
|
help="Fixed grid points where damping functions are calculated ")
|
|
|
|
parser.add_option("--gruneisen", dest="is_gruneisen",
|
|
|
|
action="store_true",
|
|
|
|
help="Calculate phonon Gruneisen parameter")
|
|
|
|
parser.add_option("--jdos",
|
|
|
|
dest="is_joint_dos",
|
|
|
|
action="store_true",
|
|
|
|
help="Calculate joint density of states")
|
2013-02-25 13:42:34 +08:00
|
|
|
parser.add_option("--lw", "--linewidth",
|
|
|
|
dest="is_linewidth",
|
2013-02-08 17:34:46 +08:00
|
|
|
action="store_true",
|
|
|
|
help="Calculate linewidths")
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--mesh",
|
|
|
|
dest="mesh_numbers",
|
|
|
|
type="string",
|
|
|
|
help="Mesh numbers")
|
2013-02-25 13:42:34 +08:00
|
|
|
parser.add_option("--md", "--mesh_divisors",
|
|
|
|
dest="mesh_divisors",
|
|
|
|
type="string",
|
|
|
|
help="Divisors for mesh numbers")
|
2013-03-12 02:27:13 +08:00
|
|
|
parser.add_option("--multiple_sigmas", dest="multiple_sigmas",
|
|
|
|
type="string",
|
|
|
|
help="Multiple sigmas for smearing width used for limited functions")
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--nac", dest="is_nac",
|
|
|
|
action="store_true",
|
|
|
|
help="Non-analytical term correction")
|
2013-02-18 14:05:15 +08:00
|
|
|
parser.add_option("--noks", "--no_kappa_stars", dest="no_kappa_stars",
|
|
|
|
action="store_true",
|
|
|
|
help="Deactivate summation of partial kappa at q-stars"),
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--nosym", dest="is_nosym",
|
|
|
|
action="store_true",
|
|
|
|
help="No symmetrization of triplets")
|
2013-02-20 17:21:35 +08:00
|
|
|
parser.add_option("--freq_step", dest="freq_step", type="float",
|
2012-12-11 16:19:23 +08:00
|
|
|
help="Pitch of frequency calculated")
|
|
|
|
parser.add_option("--pa", "--primitive_axis", dest="primitive_axis",
|
|
|
|
action="store", type="string",
|
|
|
|
help="Same as PRIMITIVE_AXIS tags")
|
|
|
|
parser.add_option("--peierls", dest="is_Peierls",
|
|
|
|
action="store_true",
|
|
|
|
help="Peierls approximation")
|
|
|
|
parser.add_option("--pm", dest="is_plusminus_displacements",
|
|
|
|
action="store_true",
|
|
|
|
help="Set plus minus displacements")
|
2013-02-05 17:00:25 +08:00
|
|
|
parser.add_option("--qpoints", dest="qpoints", type="string",
|
|
|
|
help="Calculate at specified q-points")
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--q_direction",
|
|
|
|
dest="q_direction",
|
|
|
|
type="string",
|
|
|
|
help="q-vector direction at q->0 for non-analytical term correction")
|
|
|
|
parser.add_option("--r2q", dest="r2q_TI_index", type="int",
|
|
|
|
help="Index to set translational invariance in transforming fc3 from real to reciprocal space")
|
2013-03-12 02:27:13 +08:00
|
|
|
parser.add_option("--read_gammas",
|
|
|
|
dest="read_gammas",
|
|
|
|
action="store_true",
|
|
|
|
help="Read Gammas from files")
|
2012-12-11 16:19:23 +08:00
|
|
|
parser.add_option("--sigma", dest="sigma", type="float",
|
|
|
|
help="Smearing width for DOS")
|
|
|
|
parser.add_option("--sym_fc2", dest="is_symmetrize_fc2",
|
|
|
|
action="store_true",
|
|
|
|
help="Symmetrize fc2 by index exchange")
|
|
|
|
parser.add_option("--tsym", dest="is_translational_symmetry",
|
|
|
|
action="store_true",
|
|
|
|
help="Impose translational invariance condition")
|
|
|
|
parser.add_option("--sym_fc3r", dest="is_symmetrize_fc3_r",
|
|
|
|
action="store_true",
|
|
|
|
help="Symmetrize fc3 in real space by index exchange")
|
|
|
|
parser.add_option("--sym_fc3q", dest="is_symmetrize_fc3_q",
|
|
|
|
action="store_true",
|
|
|
|
help="Symmetrize fc3 in reciprocal space by index exchange")
|
|
|
|
parser.add_option("--tmax", dest="tmax", type="string",
|
|
|
|
help="Maximum calculated temperature")
|
|
|
|
parser.add_option("--tmin", dest="tmin", type="string",
|
|
|
|
help="Minimum calculated temperature")
|
|
|
|
parser.add_option("--tstep", dest="tstep", type="string",
|
|
|
|
help="Calculated temperature step")
|
|
|
|
parser.add_option("--tolerance", dest="symprec", type="float",
|
|
|
|
help="Symmetry tolerance to search")
|
|
|
|
parser.add_option("--tp",
|
|
|
|
dest="is_read_triplets",
|
|
|
|
action="store_true",
|
|
|
|
help="Read triplets")
|
|
|
|
parser.add_option("-v", "--verbose", dest="verbose",
|
|
|
|
action="store_true",
|
|
|
|
help="Detailed run-time information is displayed")
|
2013-02-18 14:05:15 +08:00
|
|
|
parser.add_option("--loglevel", dest="log_level", type="int",
|
2012-12-11 16:19:23 +08:00
|
|
|
help="Log level")
|
2013-02-20 17:21:35 +08:00
|
|
|
parser.add_option("--ts", dest="temperatures",
|
|
|
|
type="string", help="Temperatures for damping functions")
|
2013-02-18 14:05:15 +08:00
|
|
|
parser.add_option("--wgp", "--write_grid_points", dest="write_grid_points",
|
|
|
|
action="store_true",
|
|
|
|
help="Write grid address of irreducible grid points for specified mesh numbers to ir_grid_address.yaml")
|
2013-03-02 17:52:22 +08:00
|
|
|
parser.add_option("--zheev", dest="zheev_test",
|
|
|
|
type="string", help="zheev test")
|
2012-12-11 16:19:23 +08:00
|
|
|
(options, args) = parser.parse_args()
|
|
|
|
option_list = parser.option_list
|
|
|
|
|
|
|
|
# Log level
|
|
|
|
if options.verbose:
|
|
|
|
log_level = 1
|
|
|
|
else:
|
2013-02-18 14:05:15 +08:00
|
|
|
log_level = options.log_level
|
2012-12-11 16:19:23 +08:00
|
|
|
|
|
|
|
# Title
|
|
|
|
if log_level:
|
|
|
|
print_phono3py()
|
|
|
|
|
|
|
|
# Create FORCES_THIRD
|
|
|
|
if options.forces_third_mode:
|
|
|
|
displacements = parse_disp_yaml('disp.yaml')
|
|
|
|
write_FORCES_THIRD(args, displacements)
|
|
|
|
print_end()
|
|
|
|
exit(0)
|
|
|
|
|
|
|
|
# Create DELTA_FC2_SETS
|
|
|
|
if options.delta_fc2_sets_mode:
|
|
|
|
displacements = parse_disp_yaml('disp.yaml')
|
|
|
|
write_DELTA_FC2_SETS(args, displacements)
|
|
|
|
print_end()
|
|
|
|
exit(0)
|
|
|
|
|
|
|
|
# Create FC2_FOURTH_SETS
|
|
|
|
if options.fc2_fourth_sets_mode:
|
|
|
|
displacements = parse_disp_yaml('disp.yaml')
|
|
|
|
write_FC2_FOURTH_SETS(args, displacements)
|
|
|
|
print_end()
|
|
|
|
exit(0)
|
|
|
|
|
|
|
|
# Import input files
|
|
|
|
if len(args) > 0:
|
|
|
|
phono3py_conf = Phono3pyConfParser(filename=args[0],
|
|
|
|
options=options,
|
|
|
|
option_list=option_list)
|
|
|
|
settings = phono3py_conf.get_settings()
|
|
|
|
|
|
|
|
else:
|
|
|
|
phono3py_conf = Phono3pyConfParser(options=options,
|
|
|
|
option_list=option_list)
|
|
|
|
settings = phono3py_conf.get_settings()
|
|
|
|
|
|
|
|
# Read POSCAR
|
|
|
|
if options.cell_poscar == None:
|
|
|
|
if os.path.exists('POSCAR'):
|
|
|
|
unitcell_filename = 'POSCAR'
|
|
|
|
else:
|
|
|
|
print_error("POSCAR could not be found.")
|
2013-02-18 14:05:15 +08:00
|
|
|
if log_level:
|
2012-12-11 16:19:23 +08:00
|
|
|
print_end()
|
|
|
|
sys.exit(1)
|
|
|
|
else:
|
|
|
|
if os.path.exists(options.cell_poscar):
|
|
|
|
unitcell_filename = options.cell_poscar
|
|
|
|
else:
|
|
|
|
print_error("The file \'%s\' could not be found." %
|
|
|
|
options.cell_poscar)
|
2013-02-18 14:05:15 +08:00
|
|
|
if log_level:
|
2012-12-11 16:19:23 +08:00
|
|
|
print_end()
|
|
|
|
sys.exit(1)
|
|
|
|
|
|
|
|
unitcell = read_vasp(unitcell_filename,
|
|
|
|
settings.get_chemical_symbols())
|
|
|
|
|
|
|
|
# Supercell and Symmetry
|
|
|
|
supercell = get_supercell(unitcell, settings.get_supercell_matrix())
|
|
|
|
symmetry = Symmetry(supercell, options.symprec)
|
|
|
|
|
|
|
|
if not settings.get_supercell_matrix_extra()==None:
|
|
|
|
supercell_extra = get_supercell(unitcell, settings.get_supercell_matrix_extra())
|
|
|
|
symmetry_extra = Symmetry(supercell_extra, options.symprec)
|
|
|
|
|
|
|
|
# Log
|
|
|
|
if log_level:
|
|
|
|
if options.is_translational_symmetry:
|
|
|
|
print "Translational symmetry:", options.is_translational_symmetry
|
|
|
|
if options.is_symmetrize_fc2:
|
|
|
|
print "FC2 symmetry of index exchange:", options.is_symmetrize_fc2
|
|
|
|
if options.is_symmetrize_fc3_r:
|
|
|
|
print "FC3 symmetry of index exchange in real space:", options.is_symmetrize_fc3_r
|
|
|
|
if options.is_symmetrize_fc3_q:
|
|
|
|
print "FC3 symmetry of index exchange in reciprocal space:", options.is_symmetrize_fc3_q
|
|
|
|
if not options.supercell_dimension_extra == None:
|
|
|
|
print "Extra supercell for fc2 is supplied."
|
|
|
|
if options.is_nac:
|
|
|
|
print "Non-analytical term correction:", options.is_nac
|
|
|
|
print "Spacegroup: ", symmetry.get_international_table()
|
|
|
|
|
|
|
|
###############################################################
|
|
|
|
# Create supercells with displacements and exit (pre-process) #
|
|
|
|
###############################################################
|
|
|
|
if options.is_displacement:
|
|
|
|
dds = get_third_order_displacements(
|
|
|
|
supercell,
|
|
|
|
symmetry,
|
|
|
|
is_plusminus=settings.get_is_plusminus_displacement(),
|
|
|
|
is_diagonal=settings.get_is_diagonal_displacement())
|
|
|
|
write_supercells_with_displacements(supercell,
|
|
|
|
dds,
|
|
|
|
options.amplitude)
|
|
|
|
|
|
|
|
#########################################
|
|
|
|
# Calculate third-order force constants #
|
|
|
|
#########################################
|
|
|
|
else:
|
|
|
|
primitive = Primitive(
|
|
|
|
supercell,
|
|
|
|
np.dot(np.linalg.inv(settings.get_supercell_matrix()),
|
|
|
|
settings.get_primitive_matrix()),
|
|
|
|
options.symprec)
|
|
|
|
|
|
|
|
if not settings.get_supercell_matrix_extra()==None:
|
|
|
|
primitive_extra = Primitive(
|
|
|
|
supercell_extra,
|
|
|
|
np.dot(np.linalg.inv(settings.get_supercell_matrix_extra()),
|
|
|
|
settings.get_primitive_matrix()),
|
|
|
|
options.symprec)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
if log_level:
|
|
|
|
print "------------------------ primitive cell for fc ---------------------------"
|
|
|
|
print_cell(primitive)
|
|
|
|
print "-------------------------- supercell for fc ------------------------------"
|
|
|
|
print_cell(supercell, mapping=primitive.get_supercell_to_primitive_map())
|
|
|
|
print "----------------- ratio (supercell for fc)/(primitive) -------------------"
|
|
|
|
for vec in np.dot(supercell.get_cell(), np.linalg.inv(primitive.get_cell())):
|
|
|
|
print "%5.2f"*3 % tuple(vec)
|
|
|
|
if log_level and (not settings.get_supercell_matrix_extra()==None):
|
|
|
|
print "------------------------- primitive cell extra ----------------------------"
|
|
|
|
print_cell(primitive_extra)
|
|
|
|
print "--------------------------- supercell extra -------------------------------"
|
|
|
|
print_cell(supercell_extra, mapping=primitive_extra.get_supercell_to_primitive_map())
|
|
|
|
print "--------------- ratio (supercell extra)/(primitive extra) ----------------"
|
|
|
|
for vec in np.dot(supercell_extra.get_cell(),
|
|
|
|
np.linalg.inv(primitive_extra.get_cell())):
|
|
|
|
print "%5.2f"*3 % tuple(vec)
|
|
|
|
|
2013-02-18 14:05:15 +08:00
|
|
|
# Write ir-grid points
|
|
|
|
if options.write_grid_points:
|
|
|
|
print "---------------------------------------------------------------------------"
|
|
|
|
mesh = settings.get_mesh_numbers()
|
|
|
|
if mesh is None:
|
2013-02-25 13:42:34 +08:00
|
|
|
print "To write grid points, mesh numbers have to be specified."
|
2013-02-18 14:05:15 +08:00
|
|
|
else:
|
2013-02-25 13:42:34 +08:00
|
|
|
print "Ir-grid points are written into ir_grid_points.yaml."
|
|
|
|
grid_points, weights = get_ir_grid_points(mesh, primitive)
|
|
|
|
w = open("ir_grid_points.yaml", 'w')
|
2013-02-18 14:05:15 +08:00
|
|
|
w.write("mesh: [ %d, %d, %d ]\n" % tuple(mesh))
|
2013-02-25 13:42:34 +08:00
|
|
|
if settings.get_mesh_divisors() is None:
|
|
|
|
w.write("mesh_divisors: [ 1, 1, 1 ]\n")
|
|
|
|
else:
|
|
|
|
w.write("mesh_divisors: [ %d, %d, %d ]\n" %
|
|
|
|
tuple(settings.get_mesh_divisors()))
|
|
|
|
grid_points, weights = reduce_grid_points(
|
|
|
|
settings.get_mesh_divisors(),
|
|
|
|
get_grid_address(mesh),
|
|
|
|
grid_points,
|
|
|
|
weights)
|
|
|
|
w.write("num_reduced_ir_grid_points: %d\n" % len(grid_points))
|
|
|
|
w.write("ir_grid_points: # [address, weight]\n")
|
|
|
|
for g, weight in zip(grid_points, weights):
|
2013-02-18 14:05:15 +08:00
|
|
|
w.write("- [ %d, %d ]\n" % (g, weight))
|
|
|
|
sys.exit(0)
|
|
|
|
|
2012-12-11 16:19:23 +08:00
|
|
|
# fc2
|
|
|
|
if options.read_fc2 or options.read_delta_fc2:
|
|
|
|
if log_level:
|
|
|
|
print "----- Read fc2 -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
fc2_with_dim = parse_fc2(supercell.get_number_of_atoms())
|
|
|
|
else:
|
|
|
|
if log_level:
|
|
|
|
print "----- Solve fc2 -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
forces_second = parse_FORCE_SETS(supercell.get_number_of_atoms())
|
|
|
|
fc2_with_dim = get_force_constants(forces_second,
|
|
|
|
symmetry,
|
|
|
|
supercell)
|
|
|
|
|
|
|
|
write_fc2_dat(fc2_with_dim)
|
|
|
|
|
|
|
|
if settings.get_supercell_matrix_extra()==None:
|
|
|
|
fc2 = fc2_with_dim
|
|
|
|
else:
|
|
|
|
# fc2 extra (FORCE_SETS_EXTRA)
|
|
|
|
if options.read_fc2_extra:
|
|
|
|
if log_level:
|
|
|
|
print "----- Read fc2 extra -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
fc2 = parse_fc2(supercell_extra.get_number_of_atoms(),
|
|
|
|
filename='fc2_extra.dat')
|
|
|
|
else:
|
|
|
|
if log_level:
|
|
|
|
print "----- Solve fc2 extra -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
forces_second_extra = parse_FORCE_SETS(
|
|
|
|
supercell_extra.get_number_of_atoms(),
|
|
|
|
filename="FORCE_SETS_EXTRA")
|
|
|
|
fc2 = get_force_constants(forces_second_extra,
|
|
|
|
symmetry_extra,
|
|
|
|
supercell_extra)
|
|
|
|
write_fc2_dat(fc2, filename='fc2_extra.dat')
|
|
|
|
|
|
|
|
if options.is_symmetrize_fc2:
|
|
|
|
set_permutation_symmetry(fc2)
|
|
|
|
|
|
|
|
if settings.get_is_nac():
|
|
|
|
if os.path.exists('BORN'):
|
|
|
|
if settings.get_supercell_matrix_extra()==None:
|
|
|
|
nac_params = parse_BORN(primitive)
|
|
|
|
else:
|
|
|
|
nac_params = parse_BORN(primitive_extra)
|
|
|
|
|
|
|
|
nac_q_direction = settings.get_q_direction()
|
|
|
|
else:
|
|
|
|
print_error("BORN not found")
|
2013-02-18 14:05:15 +08:00
|
|
|
if log_level:
|
2012-12-11 16:19:23 +08:00
|
|
|
print_end()
|
|
|
|
sys.exit(1)
|
|
|
|
else:
|
|
|
|
nac_params = None
|
|
|
|
nac_q_direction = None
|
|
|
|
|
|
|
|
# fc3
|
2013-03-12 02:27:13 +08:00
|
|
|
if options.is_joint_dos or options.read_gammas:
|
|
|
|
fc3 = None
|
|
|
|
else:
|
2012-12-11 16:19:23 +08:00
|
|
|
if options.read_fc3: # Read FORCES_THIRD
|
|
|
|
if log_level:
|
|
|
|
print "----- Read fc3 -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
fc3 = parse_fc3(supercell.get_number_of_atoms())
|
|
|
|
else: # fc3 from FORCES_THIRD and FORCES_SECOND
|
|
|
|
if log_level:
|
|
|
|
print "----- Solve fc3 -----"
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
|
|
|
displacements = parse_disp_yaml('disp.yaml')
|
|
|
|
if options.read_delta_fc2:
|
|
|
|
parse_DELTA_FC2_SETS(displacements)
|
|
|
|
else:
|
|
|
|
parse_DELTA_FORCES(displacements)
|
|
|
|
fc3 = get_fc3(supercell,
|
|
|
|
displacements,
|
|
|
|
fc2_with_dim,
|
|
|
|
symmetry,
|
|
|
|
options.is_translational_symmetry,
|
|
|
|
verbose=log_level)
|
|
|
|
|
|
|
|
num_atom = fc3.shape[0]
|
|
|
|
maxval1 = 0
|
|
|
|
maxval2 = 0
|
|
|
|
maxval3 = 0
|
|
|
|
for i in range(num_atom):
|
|
|
|
for j in range(num_atom):
|
|
|
|
for k in range(3):
|
|
|
|
for l in range(3):
|
|
|
|
for m in range(3):
|
|
|
|
val1 = fc3[:,i,j,k,l,m].sum()
|
|
|
|
val2 = fc3[i,:,j,k,l,m].sum()
|
|
|
|
val3 = fc3[i,j,:,k,l,m].sum()
|
|
|
|
if abs(val1) > abs(maxval1):
|
|
|
|
maxval1 = val1
|
|
|
|
if abs(val2) > abs(maxval2):
|
|
|
|
maxval2 = val2
|
|
|
|
if abs(val3) > abs(maxval3):
|
|
|
|
maxval3 = val3
|
|
|
|
|
|
|
|
print "max:", maxval1, maxval2, maxval3
|
|
|
|
|
|
|
|
|
|
|
|
# Symmetrize fc3_r
|
|
|
|
if options.is_symmetrize_fc3_r:
|
|
|
|
if log_level:
|
|
|
|
print "----- Symmetrize fc3 real space -----"
|
|
|
|
symmetrize_fc3(fc3)
|
|
|
|
|
|
|
|
# Write fc3
|
|
|
|
if not options.read_fc3:
|
|
|
|
if log_level:
|
|
|
|
print "----- Write fc3.dat -----"
|
|
|
|
write_fc3_dat(fc3)
|
|
|
|
|
|
|
|
if options.read_fc2_fourth_sets:
|
|
|
|
displacements = parse_disp_yaml('disp.yaml')
|
|
|
|
parse_FC2_FOURTH_SETS(displacements)
|
|
|
|
|
|
|
|
count = 0
|
|
|
|
for first_disp in displacements['first_atoms']:
|
|
|
|
count += 1
|
|
|
|
print count
|
|
|
|
print first_disp['fc2'][0, 0]
|
|
|
|
|
|
|
|
for first_disp in displacements['first_atoms']:
|
|
|
|
for second_disp in first_disp['second_atoms']:
|
|
|
|
for disp in second_disp['delta_fc2']:
|
|
|
|
count += 1
|
|
|
|
print count
|
|
|
|
print disp[0, 0]
|
|
|
|
|
|
|
|
|
|
|
|
sys.exit(0)
|
|
|
|
|
|
|
|
#============================
|
|
|
|
# Phonon Gruneisen parameter
|
|
|
|
#============================
|
|
|
|
if options.is_gruneisen:
|
2013-02-06 16:48:58 +08:00
|
|
|
mesh = settings.get_mesh_numbers()
|
|
|
|
|
|
|
|
if log_level:
|
|
|
|
print "------ Phonon Gruneisen parameter ------"
|
|
|
|
if mesh is not None:
|
|
|
|
print "Mesh:", mesh
|
|
|
|
sys.stdout.flush()
|
|
|
|
|
2013-02-05 17:00:25 +08:00
|
|
|
gr = get_gruneisen_parameters(fc2,
|
|
|
|
fc3,
|
|
|
|
supercell,
|
|
|
|
primitive,
|
2013-02-06 00:17:30 +08:00
|
|
|
factor=VaspToTHz,
|
|
|
|
is_ion_clamped=False,
|
2013-02-05 17:00:25 +08:00
|
|
|
symprec=options.symprec)
|
|
|
|
qpoints = settings.get_qpoints()
|
2013-02-06 16:48:58 +08:00
|
|
|
if mesh is not None:
|
|
|
|
gr.set_sampling_mesh(mesh,
|
|
|
|
is_gamma_center=True)
|
|
|
|
else:
|
|
|
|
gr.set_qpoints(qpoints)
|
2013-02-05 17:00:25 +08:00
|
|
|
gr.run()
|
2013-02-06 00:17:30 +08:00
|
|
|
gr.write_yaml()
|
2012-12-11 16:19:23 +08:00
|
|
|
|
|
|
|
#===========================
|
|
|
|
# phonon-phonon interaction
|
|
|
|
#===========================
|
2013-02-10 14:59:14 +08:00
|
|
|
else:
|
2012-12-11 16:19:23 +08:00
|
|
|
mesh = settings.get_mesh_numbers()
|
2013-02-08 22:44:52 +08:00
|
|
|
if options.grid_points is None:
|
|
|
|
grid_points = None
|
2012-12-11 16:19:23 +08:00
|
|
|
else:
|
2013-02-06 00:17:30 +08:00
|
|
|
grid_points = np.array([int(x)
|
|
|
|
for x in options.grid_points.split()])
|
2012-12-11 16:19:23 +08:00
|
|
|
|
|
|
|
if log_level:
|
|
|
|
print "------ Phonon-phonon ------"
|
2013-02-25 13:42:34 +08:00
|
|
|
print "Mesh sampling: [ %d %d %d ]" % tuple(mesh)
|
2013-02-08 22:44:52 +08:00
|
|
|
if grid_points is not None:
|
|
|
|
print "Grid points to be calculated:", grid_points
|
2012-12-11 16:19:23 +08:00
|
|
|
sys.stdout.flush()
|
|
|
|
|
2013-02-20 11:28:32 +08:00
|
|
|
if options.factor is None:
|
2012-12-11 16:19:23 +08:00
|
|
|
factor = VaspToTHz
|
|
|
|
else:
|
|
|
|
factor = options.factor
|
|
|
|
|
|
|
|
if options.is_cm:
|
|
|
|
freq_factor = THzToCm
|
|
|
|
else:
|
|
|
|
freq_factor = 1.0
|
|
|
|
|
2013-02-20 11:28:32 +08:00
|
|
|
if settings.get_sigma() is None:
|
2012-12-11 16:19:23 +08:00
|
|
|
sigma = 0.2 * freq_factor
|
|
|
|
else:
|
|
|
|
sigma = settings.get_sigma()
|
|
|
|
|
2013-02-20 11:28:32 +08:00
|
|
|
if settings.get_multiple_sigmas() is None:
|
|
|
|
multiple_sigmas = [sigma]
|
|
|
|
else:
|
|
|
|
multiple_sigmas = settings.get_multiple_sigmas()
|
|
|
|
|
|
|
|
if settings.get_omega_step() is None:
|
2013-02-20 17:21:35 +08:00
|
|
|
freq_step = 0.1 * freq_factor
|
2012-12-11 16:19:23 +08:00
|
|
|
else:
|
2013-02-20 17:21:35 +08:00
|
|
|
freq_step = settings.get_omega_step()
|
2012-12-11 16:19:23 +08:00
|
|
|
|
2013-02-20 11:28:32 +08:00
|
|
|
if options.freq_scale is None:
|
2012-12-11 16:19:23 +08:00
|
|
|
freq_scale = 1.0
|
|
|
|
else:
|
|
|
|
freq_scale = options.freq_scale
|
|
|
|
|
2013-02-08 17:34:46 +08:00
|
|
|
if settings.get_supercell_matrix_extra() is None:
|
|
|
|
supercell_dm = supercell
|
|
|
|
primitive_dm = primitive
|
|
|
|
else:
|
|
|
|
supercell_dm = supercell_extra
|
|
|
|
primitive_dm = primitive_extra
|
2012-12-11 16:19:23 +08:00
|
|
|
|
2013-02-20 17:21:35 +08:00
|
|
|
if settings.get_temperatures() is None:
|
|
|
|
temperatures = [None]
|
|
|
|
else:
|
|
|
|
temperatures = settings.get_temperatures()
|
|
|
|
|
2013-02-08 17:34:46 +08:00
|
|
|
if options.is_joint_dos:
|
|
|
|
jdos = JointDOS(supercell_dm,
|
|
|
|
primitive_dm,
|
|
|
|
mesh,
|
|
|
|
fc2,
|
|
|
|
nac_params,
|
|
|
|
sigma,
|
2013-02-20 17:21:35 +08:00
|
|
|
freq_step,
|
2013-02-08 17:34:46 +08:00
|
|
|
factor,
|
|
|
|
freq_factor,
|
|
|
|
freq_scale,
|
|
|
|
options.is_nosym,
|
|
|
|
options.symprec,
|
|
|
|
log_level)
|
2012-12-11 16:19:23 +08:00
|
|
|
jdos.get_jointDOS(grid_points,
|
2013-02-08 17:34:46 +08:00
|
|
|
filename=options.output_filename)
|
2012-12-11 16:19:23 +08:00
|
|
|
else:
|
|
|
|
phono3py = Phono3py(supercell,
|
|
|
|
primitive,
|
|
|
|
mesh,
|
|
|
|
fc3,
|
2013-02-15 00:00:13 +08:00
|
|
|
factor=factor,
|
|
|
|
freq_factor=freq_factor,
|
|
|
|
is_nosym=options.is_nosym,
|
|
|
|
is_symmetrize_fc3_q=options.is_symmetrize_fc3_q,
|
|
|
|
is_read_triplets=options.is_read_triplets,
|
|
|
|
r2q_TI_index=options.r2q_TI_index,
|
|
|
|
is_Peierls=options.is_Peierls,
|
|
|
|
symprec=options.symprec,
|
|
|
|
log_level=log_level)
|
2012-12-11 16:19:23 +08:00
|
|
|
|
2013-02-08 17:34:46 +08:00
|
|
|
phono3py.set_dynamical_matrix(fc2,
|
|
|
|
supercell_dm,
|
|
|
|
primitive_dm,
|
|
|
|
nac_params,
|
2013-02-15 00:00:13 +08:00
|
|
|
nac_q_direction,
|
|
|
|
frequency_scale_factor=freq_scale)
|
2012-12-11 16:19:23 +08:00
|
|
|
|
2013-03-02 17:52:22 +08:00
|
|
|
if options.zheev_test:
|
|
|
|
phono3py.solve_dynamical_matrix(
|
|
|
|
[float(x) for x in options.zheev_test.split()])
|
|
|
|
|
2012-12-11 16:19:23 +08:00
|
|
|
if options.is_decay_channel:
|
|
|
|
phono3py.get_decay_channels(grid_points,
|
|
|
|
settings.get_band_indices(),
|
|
|
|
options.temperature)
|
|
|
|
|
2013-02-08 17:34:46 +08:00
|
|
|
elif settings.get_is_linewidth():
|
|
|
|
phono3py.get_linewidth(
|
2013-02-06 23:29:22 +08:00
|
|
|
grid_points,
|
|
|
|
settings.get_band_indices(),
|
2013-02-20 11:28:32 +08:00
|
|
|
sigmas=multiple_sigmas,
|
2013-02-08 17:34:46 +08:00
|
|
|
t_max=settings.get_max_temperature(),
|
|
|
|
t_min=settings.get_min_temperature(),
|
|
|
|
t_step=settings.get_temperature_step(),
|
2013-02-06 23:29:22 +08:00
|
|
|
gamma_option=options.gamma_option,
|
|
|
|
filename=options.output_filename)
|
2012-12-11 16:19:23 +08:00
|
|
|
|
2013-02-18 14:05:15 +08:00
|
|
|
elif settings.get_is_bterta():
|
2013-02-10 00:30:08 +08:00
|
|
|
phono3py.get_thermal_conductivity(
|
2013-02-20 17:21:35 +08:00
|
|
|
sigmas=multiple_sigmas,
|
2013-02-08 22:44:52 +08:00
|
|
|
t_max=settings.get_max_temperature(),
|
|
|
|
t_min=settings.get_min_temperature(),
|
|
|
|
t_step=settings.get_temperature_step(),
|
2013-02-14 00:33:01 +08:00
|
|
|
grid_points=grid_points,
|
2013-02-25 13:42:34 +08:00
|
|
|
mesh_divisors=settings.get_mesh_divisors(),
|
2013-02-18 14:05:15 +08:00
|
|
|
no_kappa_stars=settings.get_no_kappa_stars(),
|
2013-03-12 02:27:13 +08:00
|
|
|
read_gammas=options.read_gammas,
|
2013-02-08 22:44:52 +08:00
|
|
|
gamma_option=options.gamma_option,
|
|
|
|
filename=options.output_filename)
|
2012-12-11 16:19:23 +08:00
|
|
|
else:
|
2013-02-07 22:09:36 +08:00
|
|
|
phono3py.get_damping_function(
|
2013-02-08 17:34:46 +08:00
|
|
|
grid_points,
|
|
|
|
settings.get_band_indices(),
|
2013-02-20 17:21:35 +08:00
|
|
|
sigmas=multiple_sigmas,
|
|
|
|
frequency_step=freq_step,
|
|
|
|
temperatures=temperatures,
|
2013-02-07 22:09:36 +08:00
|
|
|
filename=options.output_filename,
|
|
|
|
gamma_option=options.gamma_option)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
2012-12-11 16:19:23 +08:00
|
|
|
|
|
|
|
|
|
|
|
if log_level:
|
|
|
|
print_end()
|