IterMesh is now default for tdisp and tdispmat

This commit is contained in:
Atsushi Togo 2017-09-19 22:37:55 +09:00
parent 24d7c687aa
commit 5b96dfd9de
5 changed files with 92 additions and 104 deletions

View File

@ -852,37 +852,44 @@ class Phonopy(object):
are ignored.
direction:
Projection direction in reduced coordinates
Projection direction in reduced coordinates (
"""
self._thermal_displacements = None
if self._mesh is None:
print("Warning: \'set_mesh\' has to finish correctly "
"before \'set_thermal_displacements\'.")
return False
if self._mesh is not None:
eigvecs = self._mesh.get_eigenvectors()
mesh_nums = self._mesh.get_mesh_numbers()
if eigvecs is None:
print("Warning: Eigenvectors have to be calculated.")
return False
if np.prod(mesh_nums) != len(eigvecs):
print("Warning: Sampling mesh must not be symmetrized.")
return False
eigvecs = self._mesh.get_eigenvectors()
frequencies = self._mesh.get_frequencies()
mesh_nums = self._mesh.get_mesh_numbers()
iter_phonons = self._mesh
else:
if self._iter_mesh is not None:
iter_phonons = self._iter_mesh
else:
print("Warning: \'set_mesh\' has to finish correctly "
"before \'set_thermal_displacements\'.")
return False
if self._mesh.get_eigenvectors() is None:
print("Warning: Eigenvectors have to be calculated.")
return False
if direction is not None:
projection_direction = np.dot(direction, self._primitive.get_cell())
td = ThermalDisplacements(iter_phonons,
self._primitive.get_masses(),
projection_direction=projection_direction,
cutoff_frequency=cutoff_frequency)
else:
td = ThermalDisplacements(iter_phonons,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency)
if np.prod(mesh_nums) != len(eigvecs):
print("Warning: Sampling mesh must not be symmetrized.")
return False
td = ThermalDisplacements(frequencies,
eigvecs,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency)
if temperatures is None:
td.set_temperature_range(t_min, t_max, t_step)
else:
td.set_temperatures(temperatures)
if direction is not None:
td.project_eigenvectors(direction, self._primitive.get_cell())
td.run()
self._thermal_displacements = td
@ -937,23 +944,21 @@ class Phonopy(object):
print("Warning: Sampling mesh must not be symmetrized.")
return False
tdm = ThermalDisplacementMatrices(
self._mesh,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency,
lattice=self._primitive.get_cell().T)
iter_phonons = self._mesh
else:
if self._iter_mesh is not None:
tdm = ThermalDisplacementMatrices(
self._iter_mesh,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency,
lattice=self._primitive.get_cell().T)
iter_phonons = self._iter_mesh
else:
print("Warning: \'set_mesh\' has to finish correctly "
"before \'set_thermal_displacement_matrices\'.")
return False
tdm = ThermalDisplacementMatrices(
iter_phonons,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency,
lattice=self._primitive.get_cell().T)
if t_cif is None:
tdm.set_temperature_range(t_min, t_max, t_step)
else:

View File

@ -923,7 +923,6 @@ class PhonopySettings(Settings):
self._is_group_velocity = False
self._is_gamma_center = False
self._is_hdf5 = False
self._is_iter_mesh = False
self._is_little_cogroup = False
self._is_moment = False
self._is_plusminus_displacement = 'auto'
@ -1083,12 +1082,6 @@ class PhonopySettings(Settings):
def get_is_group_velocity(self):
return self._is_group_velocity
def set_is_iter_mesh(self, is_iter_mesh):
self._is_iter_mesh = is_iter_mesh
def get_is_iter_mesh(self):
return self._is_iter_mesh
def set_is_little_cogroup(self, is_little_cogroup):
self._is_little_cogroup = is_little_cogroup
@ -1329,10 +1322,6 @@ class PhonopyConfParser(ConfParser):
if opt_tdm_cif:
self._confs['tdispmat_cif'] = opt_tdm_cif
if opt.dest == 'is_iter_mesh':
if self._options.is_iter_mesh:
self._confs['iter_mesh'] = '.true.'
if opt.dest == 'projection_direction':
opt_proj_dir = self._options.projection_direction
if opt_proj_dir is not None:
@ -1575,11 +1564,6 @@ class PhonopyConfParser(ConfParser):
if len(atom_pairs) > 0:
self.set_parameter('tdistance', atom_pairs)
# Use IterMesh instead of Mesh
if conf_key == 'iter_mesh':
if confs['iter_mesh'].lower() == '.true.':
self.set_parameter('iter_mesh', True)
# Projection direction used for thermal displacements and PDOS
if conf_key == 'projection_direction':
vals = [float(x) for x in confs['projection_direction'].split()]
@ -1860,9 +1844,6 @@ class PhonopyConfParser(ConfParser):
self._settings.set_thermal_displacement_matrix_temperature(
params['tdispmat_cif'])
if 'iter_mesh' in params:
self._settings.set_is_iter_mesh(params['iter_mesh'])
# Thermal distances
if 'tdistance' in params:
self._settings.set_is_thermal_distances(True)

View File

@ -135,7 +135,6 @@ class Mesh(MeshBase):
def __next__(self):
if self._eigenvectors is None:
print("Eigenvectors was not found.")
return StopIteration
if self._q_count == len(self._qpoints):

View File

@ -40,8 +40,6 @@ from phonopy.interface.cif import write_cif_P1
class ThermalMotion(object):
def __init__(self,
frequencies, # have to be supplied in THz
eigenvectors,
masses,
cutoff_frequency=None):
@ -50,13 +48,8 @@ class ThermalMotion(object):
else:
self._cutoff_frequency = cutoff_frequency
self._distances = None
self._displacements = None
self._frequencies = frequencies
self._p_eigenvectors = None
self._eigenvectors = eigenvectors
self._masses = masses * AMU
self._masses3 = np.array([[m] * 3 for m in masses]).flatten() * AMU
self._masses3 = np.array([[m] * 3 for m in masses]).ravel() * AMU
self._temperatures = None
def get_Q2(self, freq, t): # freq in THz
@ -96,30 +89,6 @@ class ThermalMotion(object):
condition = np.logical_not(t_array < 0)
self._temperatures = np.extract(condition, t_array)
def project_eigenvectors(self, direction, lattice=None):
"""
direction
Without specifying lattice:
Projection direction in Cartesian coordinates
With lattice:
Projection direction in fractional coordinates
"""
if lattice is not None:
projector = np.dot(direction, lattice)
else:
projector = np.array(direction, dtype=float)
projector /= np.linalg.norm(projector)
self._p_eigenvectors = []
for vecs_q in self._eigenvectors:
p_vecs_q = []
for vecs in vecs_q.T:
p_vecs_q.append(np.dot(vecs.reshape(-1, 3), projector))
self._p_eigenvectors.append(np.transpose(p_vecs_q))
self._p_eigenvectors = np.array(self._p_eigenvectors)
def _get_population(self, freq, t): # freq in THz
if t < 1: # temperatue less than 1 K is approximated as 0 K.
return 0
@ -128,14 +97,29 @@ class ThermalMotion(object):
class ThermalDisplacements(ThermalMotion):
def __init__(self,
frequencies, # Have to be supplied in THz
eigenvectors,
iter_phonons,
masses,
projection_direction=None,
cutoff_frequency=None):
"""Calculate mean square displacements
iter_phonons: Mesh or IterMesh instance
masses: Atomic masses
projection_direction:
Eigenvector projection direction in Cartesian coordinates.
If None, eigenvector is not projected.
"""
self._iter_phonons = iter_phonons
if projection_direction is None:
self._projection_direction = None
else:
self._projection_direction = (projection_direction /
np.linalg.norm(projection_direction))
ThermalMotion.__init__(self,
frequencies,
eigenvectors,
masses,
cutoff_frequency=cutoff_frequency)
@ -145,24 +129,29 @@ class ThermalDisplacements(ThermalMotion):
return (self._temperatures, self._displacements)
def run(self):
freqs = self._frequencies
temps = self._temperatures
if self._p_eigenvectors is not None:
if self._projection_direction is not None:
masses = self._masses
eigvecs = self._p_eigenvectors
else:
masses = self._masses3
eigvecs = self._eigenvectors
temps = self._temperatures
disps = np.zeros((len(temps), len(masses)), dtype=float)
for fs, vecs2 in zip(freqs, abs(eigvecs) ** 2):
for f, v2 in zip(fs, vecs2.T):
for count, (fs, vecs) in enumerate(self._iter_phonons):
if self._projection_direction is not None:
p_vecs = []
for v in vecs.T:
p_vecs.append(np.dot(v.reshape(-1, 3),
self._projection_direction))
vecs2 = np.abs(p_vecs) ** 2
else:
vecs2 = (abs(vecs) ** 2).T
for f, v2 in zip(fs, vecs2):
if f > self._cutoff_frequency:
c = v2 / masses
for i, t in enumerate(temps):
disps[i] += self.get_Q2(f, t) * c
self._displacements = disps / len(freqs)
self._displacements = disps / (count + 1)
def write_yaml(self):
natom = len(self._masses)
@ -192,6 +181,19 @@ class ThermalDisplacements(ThermalMotion):
if is_legend:
pyplot.legend(plots, labels, loc='upper left')
def _project_eigenvectors(self):
"""Eigenvectors are projected along Cartesian direction"""
self._p_eigenvectors = []
for vecs_q in self._eigenvectors:
p_vecs_q = []
for vecs in vecs_q.T:
p_vecs_q.append(np.dot(vecs.reshape(-1, 3),
self._projection_direction))
self._p_eigenvectors.append(np.transpose(p_vecs_q))
self._p_eigenvectors = np.array(self._p_eigenvectors)
class ThermalDisplacementMatrices(ThermalMotion):
def __init__(self,
iter_phonons,
@ -200,8 +202,6 @@ class ThermalDisplacementMatrices(ThermalMotion):
lattice=None): # column vectors in real space
ThermalMotion.__init__(self,
None,
None,
masses,
cutoff_frequency=cutoff_frequency)
@ -297,13 +297,16 @@ class ThermalDistances(ThermalMotion):
self._supercell = supercell
self._qpoints = qpoints
self._symprec = symprec
self._frequencies = frequencies
self._eigenvectors = eigenvectors
ThermalMotion.__init__(self,
frequencies,
eigenvectors,
primitive.get_masses(),
cutoff_frequency=cutoff_frequency)
self._p_eigenvectors = None
self._distances = None
def _get_cross(self, v, delta_r, q, atom1, atom2):
phase = np.exp(2j * np.pi * np.dot(delta_r, q))
cross_val = v[atom1]*phase*v[atom2].conjugate()

View File

@ -732,8 +732,8 @@ elif run_mode == 'band' or run_mode == 'mesh' or run_mode == 'band_mesh':
q_symmetry,
is_gamma_center) = settings.get_mesh()
if settings.get_is_iter_mesh():
print("Use IterMesh")
if (settings.get_is_thermal_displacements() or
settings.get_is_thermal_displacement_matrices()):
phonon.set_iter_mesh(mesh,
mesh_shift,
is_time_reversal=t_symmetry,