Allow calculation of elements of kappa with specific band indices

This commit is contained in:
Atsushi Togo 2015-10-04 12:58:12 +09:00
parent 25b4f0e1a7
commit 5a86da02e2
5 changed files with 66 additions and 40 deletions

View File

@ -833,19 +833,15 @@ def write_kappa_to_hdf5(temperature,
weight=None,
mesh_divisors=None,
grid_point=None,
band_indices=None,
sigma=None,
filename=None):
suffix = "-m%d%d%d" % tuple(mesh)
if mesh_divisors is not None:
if (np.array(mesh_divisors, dtype=int) != 1).any():
suffix += "-d%d%d%d" % tuple(mesh_divisors)
if grid_point is not None:
suffix += ("-g%d" % grid_point)
if sigma is not None:
sigma_str = ("%f" % sigma).rstrip('0').rstrip('\.')
suffix += "-s" + sigma_str
if filename is not None:
suffix += "." + filename
suffix = _get_filename_suffix(mesh,
mesh_divisors=mesh_divisors,
grid_point=grid_point,
band_indices=band_indices,
sigma=sigma,
filename=filename)
w = h5py.File("kappa" + suffix + ".hdf5", 'w')
w.create_dataset('temperature', data=temperature)
if frequency is not None:
@ -1424,12 +1420,21 @@ def parse_grid_address(filename):
return np.array(grid_address)
def _get_filename_suffix(mesh,
mesh_divisors=None,
grid_point=None,
band_indices=None,
sigma=None,
filename=None):
suffix = "-m%d%d%d" % tuple(mesh)
if mesh_divisors is not None:
if (np.array(mesh_divisors, dtype=int) != 1).any():
suffix += "-d%d%d%d" % tuple(mesh_divisors)
if grid_point is not None:
suffix += ("-g%d" % grid_point)
if band_indices is not None:
suffix += "-"
for bi in band_indices:
suffix += "b%d" % (bi + 1)
if sigma is not None:
sigma_str = ("%f" % sigma).rstrip('0').rstrip('\.')
suffix += "-s" + sigma_str

View File

@ -107,6 +107,9 @@ class Phono3py:
self._band_indices_flatten = None
if mesh is not None:
self._mesh = np.array(mesh, dtype='intc')
self.set_band_indices(band_indices)
def set_band_indices(self, band_indices):
if band_indices is None:
num_band = self._primitive.get_number_of_atoms() * 3
self._band_indices = [np.arange(num_band, dtype='intc')]

View File

@ -291,7 +291,8 @@ class Conductivity:
def _set_gv(self, i):
# Group velocity [num_freqs, 3]
self._gv[i] = self._get_gv(self._qpoints[i])
gv = self._get_gv(self._qpoints[i])
self._gv[i] = gv[self._pp.get_band_indices(), :]
def _get_gv(self, q):
return get_group_velocity(

View File

@ -60,12 +60,12 @@ def get_thermal_conductivity_RTA(
if log_level > 1 and read_gamma is False:
_write_triplets(interaction)
if grid_points is None:
if (grid_points is None and _all_bands_exist(interaction)):
br.set_kappa_at_sigmas()
_write_kappa(br, filename=output_filename, log_level=log_level)
return br
def _write_gamma(br, interaction, i, filename=None):
grid_points = br.get_grid_points()
group_velocities = br.get_group_velocities()
@ -77,10 +77,15 @@ def _write_gamma(br, interaction, i, filename=None):
gamma = br.get_gamma()
gamma_isotope = br.get_gamma_isotope()
sigmas = br.get_sigmas()
gp = grid_points[i]
frequencies = interaction.get_phonons()[0][gp]
if _all_bands_exist(interaction):
band_indices = None
frequencies = interaction.get_phonons()[0][gp]
else:
band_indices = interaction.get_band_indices()
frequencies = interaction.get_phonons()[0][gp][band_indices]
for j, sigma in enumerate(sigmas):
if gamma_isotope is not None:
gamma_isotope_at_sigma = gamma_isotope[j, i]
@ -97,9 +102,18 @@ def _write_gamma(br, interaction, i, filename=None):
averaged_pp_interaction=ave_pp[i],
mesh_divisors=mesh_divisors,
grid_point=gp,
band_indices=band_indices,
sigma=sigma,
filename=filename)
def _all_bands_exist(interaction):
band_indices = interaction.get_band_indices()
num_band = interaction.get_primitive().get_number_of_atoms() * 3
if len(band_indices) == num_band:
if (band_indices - np.arange(num_band) == 0).all():
return True
return False
def _write_triplets(interaction, filename=None):
triplets, weights = interaction.get_triplets_at_q()[:2]
grid_address = interaction.get_grid_address()
@ -368,6 +382,7 @@ class Conductivity_RTA(Conductivity):
i = self._grid_point_count
self._show_log_header(i)
grid_point = self._grid_points[i]
if self._read_gamma:
if self._use_ave_pp:
self._collision.set_grid_point(grid_point)
@ -389,13 +404,15 @@ class Conductivity_RTA(Conductivity):
if self._isotope is not None and not self._read_gamma_iso:
self._set_gamma_isotope_at_sigmas(i)
self._cv[:, i, :] = self._get_cv(self._frequencies[grid_point])
freqs = self._frequencies[grid_point][self._pp.get_band_indices()]
self._cv[:, i, :] = self._get_cv(freqs)
self._set_gv(i)
if self._log_level:
self._show_log(self._qpoints[i], i)
def _allocate_values(self):
num_band0 = len(self._pp.get_band_indices())
num_band = self._primitive.get_number_of_atoms() * 3
num_grid_points = len(self._grid_points)
self._kappa = np.zeros((len(self._sigmas),
@ -404,28 +421,23 @@ class Conductivity_RTA(Conductivity):
self._mode_kappa = np.zeros((len(self._sigmas),
len(self._temperatures),
num_grid_points,
num_band,
num_band0,
6), dtype='double')
if not self._read_gamma:
self._gamma = np.zeros((len(self._sigmas),
len(self._temperatures),
num_grid_points,
num_band), dtype='double')
self._gv = np.zeros((num_grid_points,
num_band,
3), dtype='double')
self._cv = np.zeros((len(self._temperatures),
num_grid_points,
num_band), dtype='double')
num_band0), dtype='double')
self._gv = np.zeros((num_grid_points, num_band0, 3), dtype='double')
self._cv = np.zeros(
(len(self._temperatures), num_grid_points, num_band0), dtype='double')
if self._isotope is not None:
self._gamma_iso = np.zeros((len(self._sigmas),
num_grid_points,
num_band), dtype='double')
self._averaged_pp_interaction = np.zeros((num_grid_points, num_band),
dtype='double')
self._num_ignored_phonon_modes = np.zeros((len(self._sigmas),
len(self._temperatures)),
dtype='intc')
self._gamma_iso = np.zeros(
(len(self._sigmas), num_grid_points, num_band0), dtype='double')
self._averaged_pp_interaction = np.zeros(
(num_grid_points, num_band0), dtype='double')
self._num_ignored_phonon_modes = np.zeros(
(len(self._sigmas), len(self._temperatures)), dtype='intc')
self._collision = ImagSelfEnergy(
self._pp,
unit_conversion=self._gamma_unit_conversion)
@ -483,7 +495,7 @@ class Conductivity_RTA(Conductivity):
def _show_log(self, q, i):
gp = self._grid_points[i]
frequencies = self._frequencies[gp]
frequencies = self._frequencies[gp][self._pp.get_band_indices()]
gv = self._gv[i]
ave_pp = self._averaged_pp_interaction[i]

View File

@ -26,11 +26,9 @@ class Interaction:
self._primitive = primitive
self._mesh = np.array(mesh, dtype='intc')
self._symmetry = symmetry
num_band = primitive.get_number_of_atoms() * 3
if band_indices is None:
self._band_indices = np.arange(num_band, dtype='intc')
else:
self._band_indices = np.array(band_indices, dtype='intc')
self._band_indices = None
self.set_band_indices(band_indices)
self._constant_averaged_interaction = constant_averaged_interaction
self._frequency_factor_to_THz = frequency_factor_to_THz
@ -119,6 +117,13 @@ class Interaction:
def get_band_indices(self):
return self._band_indices
def set_band_indices(self, band_indices):
num_band = self._primitive.get_number_of_atoms() * 3
if band_indices is None:
self._band_indices = np.arange(num_band, dtype='intc')
else:
self._band_indices = np.array(band_indices, dtype='intc')
def get_frequency_factor_to_THz(self):
return self._frequency_factor_to_THz