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

A boundary treatment for no-kappa-stars option in thermal conductivity

This commit is contained in:
Atsushi Togo 2013-09-27 16:15:33 +09:00
parent dab35c41c9
commit 1b7e3aa0e8
8 changed files with 115 additions and 82 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

61
anharmonic/phonon3/conductivity_RTA.py
View File

@ -59,6 +59,7 @@ class conductivity_RTA:
self._grid_points = None
self._grid_weights = None
self._grid_address = None
self._grid_address_with_boundary = None # Used only when no_kappa_stars
self._gamma = None
self._read_gamma = False
@ -106,6 +107,8 @@ class conductivity_RTA:
coarse_grid_address = get_grid_address(self._coarse_mesh)
coarse_grid_points = np.arange(np.prod(self._coarse_mesh),
dtype='intc')
self._grid_address_with_boundary = get_bz_grid_address(
self._mesh, primitive_lattice, with_boundary=True)
self._grid_points = from_coarse_to_dense_grid_points(
self._mesh,
self._mesh_divisors,
@ -261,22 +264,46 @@ class conductivity_RTA:
def _get_gv_by_gv(self, i):
grid_address = self._grid_address[self._grid_points[i]]
# Group velocity [num_freqs, 3]
gv = get_group_velocity(
self._qpoint,
self._dynamical_matrix,
q_length=self._gv_delta_q,
symmetry=self._symmetry,
frequency_factor_to_THz=self._frequency_factor_to_THz)
self._gv[i] = gv
rotation_map = self._get_rotation_map_for_star(grid_address)
gv_by_gv = self._get_gv_by_gv_on_star(gv, rotation_map)
if self._no_kappa_stars:
gv_by_gv_tmp = []
rotation_map = [0]
for address in self._grid_address_with_boundary:
if ((grid_address - address) % self._mesh == 0).all():
qpoint = address.astype('double') / self._mesh
gv = get_group_velocity(
qpoint,
self._dynamical_matrix,
q_length=self._gv_delta_q,
symmetry=self._symmetry,
frequency_factor_to_THz=self._frequency_factor_to_THz)
self._gv[i] = gv
gv_by_gv_tmp.append(
self._get_gv_by_gv_on_star(gv, rotation_map)[0])
if self._log_level:
self._show_log(self._qpoint,
self._frequencies[i],
gv,
rotation_map)
if self._log_level:
self._show_log(qpoint,
self._frequencies[i],
gv,
rotation_map)
gv_by_gv = [np.sum(gv_by_gv_tmp, axis=0) / len(gv_by_gv_tmp)]
else:
# Group velocity [num_freqs, 3]
gv = get_group_velocity(
self._qpoint,
self._dynamical_matrix,
q_length=self._gv_delta_q,
symmetry=self._symmetry,
frequency_factor_to_THz=self._frequency_factor_to_THz)
self._gv[i] = gv
rotation_map = self._get_rotation_map_for_star(grid_address)
gv_by_gv = self._get_gv_by_gv_on_star(gv, rotation_map)
if self._log_level:
self._show_log(self._qpoint,
self._frequencies[i],
gv,
rotation_map)
# check if the number of rotations is correct.
if self._grid_weights is not None:
@ -419,7 +446,7 @@ class conductivity_RTA:
frequencies,
group_velocity,
rotation_map):
print "Frequency, projected group velocity (x, y, z), norm at k-stars",
print "Frequency, projected group velocity (x, y, z), group velocity norm",
if self._gv_delta_q is None:
print
else:
@ -474,5 +501,3 @@ class conductivity_RTA:
grid_address,
grid_point=grid_point,
filename=self._filename)

7
anharmonic/phonon3/triplets.py
View File

@ -53,12 +53,15 @@ def get_grid_address(mesh):
return grid_address
def get_bz_grid_address(mesh, primitive_lattice):
def get_bz_grid_address(mesh, primitive_lattice, with_boundary=False):
grid_address = get_grid_address(mesh)
bz_grid_address, bz_map = spg.relocate_BZ_grid_address(grid_address,
mesh,
primitive_lattice)
return grid_address
if with_boundary:
return bz_grid_address
else:
return grid_address
def get_grid_point_from_address(address, mesh, with_boundary=False):
# X runs first in XYZ

15
c/_spglib.c
View File

@ -554,15 +554,16 @@ static PyObject * relocate_BZ_grid_address(PyObject *self, PyObject *args)
const int* is_shift = (int*)is_shift_py->data;
SPGCONST double (*reciprocal_lattice)[3] =
(double(*)[3])reciprocal_lattice_py->data;
int num_ir_gp;
spg_relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
reciprocal_lattice,
is_shift);
num_ir_gp = spg_relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
reciprocal_lattice,
is_shift);
Py_RETURN_NONE;
return PyInt_FromLong((long) num_ir_gp);
}
static PyObject * get_triplets_reciprocal_mesh_at_q(PyObject *self, PyObject *args)

50
c/spglib/kpoint.c
View File

@ -65,12 +65,12 @@ get_ir_reciprocal_mesh_openmp(int grid_address[][3],
const int mesh[3],
const int is_shift[3],
SPGCONST PointSymmetry * point_symmetry);
static void relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
static int relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
static double get_tolerance_for_BZ_reduction(SPGCONST double rec_lattice[3][3]);
static int get_ir_triplets_at_q(int weights[],
int grid_address[][3],
@ -204,19 +204,19 @@ int kpt_get_stabilized_reciprocal_mesh(int grid_address[][3],
}
void kpt_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
int kpt_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
{
relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
rec_lattice,
is_shift);
return relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
rec_lattice,
is_shift);
}
int kpt_get_ir_triplets_at_q(int weights[],
@ -558,12 +558,12 @@ get_ir_reciprocal_mesh_openmp(int grid_address[][3],
/* Relocate grid addresses to first Brillouin zone */
/* bz_grid_address[prod(mesh + 1)][3] */
/* bz_map[prod(mesh * 2)] */
static void relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
static int relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
{
double tolerance, min_distance;
double vector[3], distance[27];
@ -618,6 +618,8 @@ static void relocate_BZ_grid_address(int bz_grid_address[][3],
grid_address[i][j] += search_space[min_index][j] * mesh[j];
}
}
return num_gp;
}
static double get_tolerance_for_BZ_reduction(SPGCONST double rec_lattice[3][3])

24
c/spglib/spglib.c
View File

@ -524,19 +524,19 @@ int spg_get_stabilized_reciprocal_mesh(int grid_address[][3],
qpoints);
}
void spg_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
int spg_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3])
{
kpt_relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
rec_lattice,
is_shift);
return kpt_relocate_BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
mesh,
rec_lattice,
is_shift);
}
int spg_get_triplets_reciprocal_mesh_at_q(int weights[],

12
c/spglib_h/kpoint.h
View File

@ -27,12 +27,12 @@ int kpt_get_stabilized_reciprocal_mesh(int grid_address[][3],
const MatINT * pointgroup_real,
const int num_q,
SPGCONST double qpoints[][3]);
void kpt_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
int kpt_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
int kpt_get_ir_triplets_at_q(int weights[],
int grid_address[][3],
int third_q[],

13
c/spglib_h/spglib.h
View File

@ -301,12 +301,13 @@ int spg_get_stabilized_reciprocal_mesh(int grid_address[][3],
SPGCONST double qpoints[][3]);
/* Grid addresses are relocated inside Brillouin zone. */
void spg_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
/* Number of ir-grid-points inside Brillouin zone is returned. */
int spg_relocate_BZ_grid_address(int bz_grid_address[][3],
int bz_map[],
int grid_address[][3],
const int mesh[3],
SPGCONST double rec_lattice[3][3],
const int is_shift[3]);
/* Irreducible triplets of k-points are searched under conservation of */
/* :math:``\mathbf{k}_1 + \mathbf{k}_2 + \mathbf{k}_3 = \mathbf{G}``. */

15
phonopy/structure/spglib.py
View File

@ -252,14 +252,15 @@ def relocate_BZ_grid_address(grid_address,
bz_grid_address = np.zeros(
((mesh[0] + 1) * (mesh[1] + 1) * (mesh[2] + 1), 3), dtype='intc')
bz_map = np.zeros(np.prod(mesh) * 8, dtype='intc')
spg.BZ_grid_address(bz_grid_address,
bz_map,
grid_address,
np.array(mesh, dtype='intc').copy(),
np.array(reciprocal_lattice, dtype='double').copy(),
np.array(is_shift, dtype='intc').copy())
num_bz_ir = spg.BZ_grid_address(
bz_grid_address,
bz_map,
grid_address,
np.array(mesh, dtype='intc').copy(),
np.array(reciprocal_lattice, dtype='double').copy(),
np.array(is_shift, dtype='intc').copy())
return bz_grid_address, bz_map
return bz_grid_address[:num_bz_ir], bz_map
def get_stabilized_reciprocal_mesh(mesh,
rotations,