scnc
/
phonopy
mirror of https://github.com/phonopy/phonopy.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Made temperatures keyword

This commit is contained in:
Atsushi Togo 2016-10-21 17:07:40 +09:00
parent 5875dc64f4
commit 51391bb192
3 changed files with 136 additions and 67 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

83
phonopy/api_phonopy.py
View File

@ -567,41 +567,6 @@ class Phonopy(object):
return plt
# Thermal property
def set_thermal_properties(self,
t_step=10,
t_max=1000,
t_min=0,
is_projection=False,
band_indices=None,
cutoff_frequency=None):
if self._mesh is None:
print("Warning: set_mesh has to be done before "
"set_thermal_properties")
return False
else:
tp = ThermalProperties(self._mesh.get_frequencies(),
weights=self._mesh.get_weights(),
eigenvectors=self._mesh.get_eigenvectors(),
is_projection=is_projection,
band_indices=band_indices,
cutoff_frequency=cutoff_frequency)
tp.run(t_step=t_step, t_max=t_max, t_min=t_min)
self._thermal_properties = tp
def get_thermal_properties(self):
temps, fe, entropy, cv = \
self._thermal_properties.get_thermal_properties()
return temps, fe, entropy, cv
def plot_thermal_properties(self):
import matplotlib.pyplot as plt
self._thermal_properties.plot(plt)
return plt
def write_yaml_thermal_properties(self, filename='thermal_properties.yaml'):
self._thermal_properties.write_yaml(filename=filename)
# DOS
def set_total_DOS(self,
sigma=None,
@ -712,11 +677,54 @@ class Phonopy(object):
def write_partial_DOS(self):
self._pdos.write()
# Thermal property
def set_thermal_properties(self,
t_step=10,
t_max=1000,
t_min=0,
temperatures=None,
is_projection=False,
band_indices=None,
cutoff_frequency=None):
if self._mesh is None:
print("Warning: set_mesh has to be done before "
"set_thermal_properties")
return False
else:
tp = ThermalProperties(self._mesh.get_frequencies(),
weights=self._mesh.get_weights(),
eigenvectors=self._mesh.get_eigenvectors(),
is_projection=is_projection,
band_indices=band_indices,
cutoff_frequency=cutoff_frequency)
if temperatures is None:
tp.set_temperature_range(t_step=t_step,
t_max=t_max,
t_min=t_min)
else:
tp.set_temperatures(temperatures)
tp.run()
self._thermal_properties = tp
def get_thermal_properties(self):
temps, fe, entropy, cv = \
self._thermal_properties.get_thermal_properties()
return temps, fe, entropy, cv
def plot_thermal_properties(self):
import matplotlib.pyplot as plt
self._thermal_properties.plot(plt)
return plt
def write_yaml_thermal_properties(self, filename='thermal_properties.yaml'):
self._thermal_properties.write_yaml(filename=filename)
# Thermal displacement
def set_thermal_displacements(self,
t_step=10,
t_max=1000,
t_min=0,
temperatures=None,
direction=None,
cutoff_frequency=None):
"""
@ -750,7 +758,10 @@ class Phonopy(object):
eigvecs,
self._primitive.get_masses(),
cutoff_frequency=cutoff_frequency)
td.set_temperature_range(t_min, t_max, t_step)
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()

54
phonopy/phonon/thermal_displacement.py
View File

@ -59,30 +59,42 @@ class ThermalMotion(object):
self._masses3 = np.array([[m] * 3 for m in masses]).flatten() * AMU
self._temperatures = None
def _get_population(self, freq, t): # freq in THz
if t < 1: # temperatue less than 1 K is approximated as 0 K.
return 0
else:
return 1.0 / (np.exp(freq * THzToEv / (Kb * t)) - 1)
def get_Q2(self, freq, t): # freq in THz
return Hbar * EV / Angstrom ** 2 * (
(self._get_population(freq, t) + 0.5) / (freq * 1e12 * 2 * np.pi))
def set_temperature_range(self, t_min=0, t_max=1000, t_step=10):
if t_min < 0:
t_min = 0
if t_step < 0:
t_step = 0
temps = []
t = t_min
while t < t_max + t_step / 2.0:
temps.append(t)
t += t_step
self._temperatures = np.array(temps)
def get_temperatures(self):
return self._temperatures
def set_temperature_range(self, t_min=None, t_max=None, t_step=None):
if t_min is None:
_t_min = 10
elif t_min < 0:
_t_min = 0
else:
_t_min = t_min
if t_max is None:
_t_max = 1000
elif t_max > _t_min:
_t_max = t_max
else:
_t_max = _t_min
if t_step is None:
_t_step = 10
elif t_step > 0:
_t_step = t_step
else:
_t_step = 10
self._temperatures = np.arange(_t_min, _t_max + _t_step / 2.0, _t_step,
dtype='double')
def set_temperatures(self, temperatures):
self._temperatures = temperatures
t_array = np.array(temperatures)
condition = np.logical_not(t_array < 0)
self._temperatures = np.extract(condition, t_array)
def project_eigenvectors(self, direction, lattice=None):
"""
@ -108,6 +120,12 @@ class ThermalMotion(object):
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
else:
return 1.0 / (np.exp(freq * THzToEv / (Kb * t)) - 1)
class ThermalDisplacements(ThermalMotion):
def __init__(self,
frequencies, # Have to be supplied in THz

66
phonopy/phonon/thermal_properties.py
View File

@ -68,7 +68,7 @@ class ThermalPropertiesBase(object):
self._weights = None
self._is_projection = is_projection
self._cutoff_frequency = cutoff_frequency
if band_indices is not None:
bi = np.hstack(band_indices).astype('intc')
self._band_indices = bi
@ -149,8 +149,17 @@ class ThermalProperties(ThermalPropertiesBase):
is_projection=is_projection,
band_indices=band_indices,
cutoff_frequency=cutoff_frequency)
self._thermal_properties = None
self._temperatures = None
self._high_T_entropy = None
self._zero_point_energy = None
self._projected_thermal_properties = None
self._set_high_T_entropy_and_zero_point_energy()
def get_temperatures(self):
return self._temperatures
def get_number_of_integrated_modes(self):
"""Number of phonon modes used for integration on sampling mesh"""
return self._num_integrated_modes
@ -165,6 +174,36 @@ class ThermalProperties(ThermalPropertiesBase):
def get_high_T_entropy(self):
return self._high_T_entropy
def set_temperature_range(self, t_min=None, t_max=None, t_step=None):
if t_min is None:
_t_min = 10
elif t_min < 0:
_t_min = 0
else:
_t_min = t_min
if t_max is None:
_t_max = 1000
elif t_max > _t_min:
_t_max = t_max
else:
_t_max = _t_min
if t_step is None:
_t_step = 10
elif t_step > 0:
_t_step = t_step
else:
_t_step = 10
self._temperatures = np.arange(_t_min, _t_max + _t_step / 2.0, _t_step,
dtype='double')
def set_temperatures(self, temperatures):
t_array = np.array(temperatures)
condition = np.logical_not(t_array < 0)
self._temperatures = np.extract(condition, t_array)
def plot(self, pyplot):
temps, fe, entropy, cv = self._thermal_properties
@ -177,34 +216,35 @@ class ThermalProperties(ThermalPropertiesBase):
pyplot.grid(True)
pyplot.xlabel('Temperature [K]')
def set_thermal_properties(self, t_step=10, t_max=1000, t_min=0):
def run(self, t_step=None, t_max=None, t_min=None):
import warnings
warnings.warn("\'set_thermal_properties\' method is depreciated. "
"Use \'run\' method instead.")
self.run(t_step=t_step, t_max=t_max, t_min=t_min)
if (t_step is not None or
t_max is not None or
t_min is not None):
warnings.warn("keywords for this method are depreciated. "
"Use \'set_temperature_range\' or "
"\'set_temperature_range\' method instead.")
self.set_temperature_range(t_min=t_min, t_max=t_max, t_step=t_step)
def run(self, t_step=10, t_max=1000, t_min=0):
temperatures = np.arange(t_min, t_max + t_step / 2.0, t_step,
dtype='double')
fe = []
entropy = []
cv = []
energy = []
try:
import phonopy._phonopy as phonoc
for t in temperatures:
for t in self._temperatures:
props = self._get_c_thermal_properties(t)
fe.append(props[0] * EvTokJmol + self._zero_point_energy)
entropy.append(props[1] * EvTokJmol * 1000)
cv.append(props[2] * EvTokJmol * 1000)
except ImportError:
for t in temperatures:
for t in self._temperatures:
props = self._get_py_thermal_properties(t)
fe.append(props[0])
entropy.append(props[1] * 1000,)
cv.append(props[2] * 1000)
self._thermal_properties = [temperatures,
self._thermal_properties = [self._temperatures,
np.array(fe, dtype='double', order='C'),
np.array(entropy, dtype='double', order='C'),
np.array(cv, dtype='double', order='C')]
@ -214,13 +254,13 @@ class ThermalProperties(ThermalPropertiesBase):
entropy = []
cv = []
energy = []
for t in temperatures:
for t in self._temperatures:
fe.append(self.get_free_energy(t))
entropy.append(self.get_entropy(t) * 1000,)
cv.append(self.get_heat_capacity_v(t) * 1000)
self._projected_thermal_properties = [
temperatures,
self._temperatures,
np.array(fe, dtype='double'),
np.array(entropy, dtype='double'),
np.array(cv, dtype='double')]