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

Refactoring of bzgrid.c

This commit is contained in:
Atsushi Togo 2021-04-24 16:14:44 +09:00
parent 5c689564b4
commit 52912e971a
6 changed files with 85 additions and 146 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

144
c/bzgrid.c
View File

@ -168,9 +168,7 @@ static long bz_search_space[BZG_NUM_BZ_SEARCH_SPACE][3] = {
{-1, -1, -1}
};
static RotMats *get_point_group_reciprocal(const RotMats * rotations,
const long is_time_reversal,
const long is_transose);
static long get_ir_grid_map(long ir_mapping_table[],
const long D_diag[3],
const long PS[3],
@ -197,24 +195,31 @@ static void multiply_matrix_vector_d3(double v[3],
const double b[3]);
static double norm_squared_d3(const double a[3]);
long bzg_get_ir_grid_map(long ir_mapping_table[],
const long D_diag[3],
const long PS[3],
const RotMats *rot_reciprocal)
{
long num_ir;
num_ir = get_ir_grid_map(ir_mapping_table,
D_diag,
PS,
rot_reciprocal);
return num_ir;
}
RotMats *bzg_get_point_group_reciprocal(const RotMats * rotations,
const long is_time_reversal)
RotMats *bzg_get_reciprocal_point_group(const long (*rec_rotations_in)[3][3],
const long num_rot,
const long is_time_reversal,
const long is_transpose)
{
return get_point_group_reciprocal(rotations, is_time_reversal, 0);
long i, num_rot_out;
long rec_rotations_out[48][3][3];
RotMats *rec_rotations;
num_rot_out = grg_get_reciprocal_point_group(rec_rotations_out,
rec_rotations_in,
num_rot,
is_time_reversal,
is_transpose);
if (num_rot_out == 0) {
return NULL;
}
rec_rotations = bzg_alloc_RotMats(num_rot_out);
for (i = 0; i < num_rot_out; i++) {
lagmat_copy_matrix_l3(rec_rotations->mat[i], rec_rotations_out[i]);
}
return rec_rotations;
}
long bzg_get_ir_reciprocal_mesh(long *ir_mapping_table,
@ -224,16 +229,17 @@ long bzg_get_ir_reciprocal_mesh(long *ir_mapping_table,
const long (*rec_rotations_in)[3][3],
const long num_rot)
{
long i, num_ir;
RotMats *rotations, *rot_reciprocal;
long num_ir;
RotMats *rot_reciprocal;
rotations = bzg_alloc_RotMats(num_rot);
for (i = 0; i < num_rot; i++) {
lagmat_copy_matrix_l3(rotations->mat[i], rec_rotations_in[i]);
rot_reciprocal = bzg_get_reciprocal_point_group(rec_rotations_in,
num_rot,
is_time_reversal,
0);
if (rot_reciprocal == NULL) {
return 0;
}
rot_reciprocal = NULL;
rot_reciprocal = get_point_group_reciprocal(rotations, is_time_reversal, 0);
num_ir = get_ir_grid_map(ir_mapping_table,
D_diag,
PS,
@ -241,8 +247,6 @@ long bzg_get_ir_reciprocal_mesh(long *ir_mapping_table,
bzg_free_RotMats(rot_reciprocal);
rot_reciprocal = NULL;
bzg_free_RotMats(rotations);
rotations = NULL;
return num_ir;
}
@ -376,90 +380,6 @@ long bzg_inverse_unimodular_matrix_l3(long m[3][3],
return det;
}
/* Return NULL if failed */
static RotMats *get_point_group_reciprocal(const RotMats * rotations,
const long is_time_reversal,
const long is_transpose)
{
long i, j, num_rot;
RotMats *rot_reciprocal, *rot_return;
long *unique_rot;
const long inversion[3][3] = {
{-1, 0, 0 },
{ 0,-1, 0 },
{ 0, 0,-1 }
};
rot_reciprocal = NULL;
rot_return = NULL;
unique_rot = NULL;
if (is_time_reversal) {
if ((rot_reciprocal = bzg_alloc_RotMats(rotations->size * 2)) == NULL) {
return NULL;
}
} else {
if ((rot_reciprocal = bzg_alloc_RotMats(rotations->size)) == NULL) {
return NULL;
}
}
if ((unique_rot = (long*)malloc(sizeof(long) * rot_reciprocal->size)) == NULL) {
warning_print("Memory of unique_rot could not be allocated.");
bzg_free_RotMats(rot_reciprocal);
rot_reciprocal = NULL;
return NULL;
}
for (i = 0; i < rot_reciprocal->size; i++) {
unique_rot[i] = -1;
}
if (is_transpose) {
for (i = 0; i < rotations->size; i++) {
lagmat_transpose_matrix_l3(rot_reciprocal->mat[i], rotations->mat[i]);
}
} else {
for (i = 0; i < rotations->size; i++) {
lagmat_copy_matrix_l3(rot_reciprocal->mat[i], rotations->mat[i]);
}
}
if (is_time_reversal) {
for (i = 0; i < rotations->size; i++) {
lagmat_multiply_matrix_l3(rot_reciprocal->mat[rotations->size + i],
inversion,
rot_reciprocal->mat[i]);
}
}
num_rot = 0;
for (i = 0; i < rot_reciprocal->size; i++) {
for (j = 0; j < num_rot; j++) {
if (lagmat_check_identity_matrix_l3(rot_reciprocal->mat[unique_rot[j]],
rot_reciprocal->mat[i])) {
goto escape;
}
}
unique_rot[num_rot] = i;
num_rot++;
escape:
;
}
if ((rot_return = bzg_alloc_RotMats(num_rot)) != NULL) {
for (i = 0; i < num_rot; i++) {
lagmat_copy_matrix_l3(rot_return->mat[i], rot_reciprocal->mat[unique_rot[i]]);
}
}
free(unique_rot);
unique_rot = NULL;
bzg_free_RotMats(rot_reciprocal);
rot_reciprocal = NULL;
return rot_return;
}
/* It is assumed that the rotations have been examined by
* grg_transform_rotations, i.e., no broken symmetry of grid is ensured. */

10
c/bzgrid.h
View File

@ -92,12 +92,10 @@ typedef struct {
} ConstBZGrid;
long bzg_get_ir_grid_map(long ir_mapping_table[],
const long D_diag[3],
const long PS[3],
const RotMats *rot_reciprocal);
RotMats *bzg_get_point_group_reciprocal(const RotMats * rotations,
const long is_time_reversal);
RotMats *bzg_get_reciprocal_point_group(const long (*rec_rotations_in)[3][3],
const long num_rot,
const long is_time_reversal,
const long is_transpose);
long bzg_get_ir_reciprocal_mesh(long *ir_mapping_table,
const long D_diag[3],
const long PS[3],

27
c/grgrid.c
View File

@ -250,14 +250,20 @@ void grg_get_ir_grid_map(long *ir_grid_map,
PS);
}
/* Extract unique rotation matrices and transpose them. */
/* When is_time_reversal == 1, inverse of the extracted matrices are */
/* included. */
/* Unique reciprocal rotations are collected from input rotations. */
/* is_transpose == 0 : Input rotations are considered those for */
/* reciprocal space. */
/* is_transpose != 0 : Input rotations are considered those for */
/* direct space, i.e., the rotation matrices are transposed. */
/* is_time_reversal controls if inversion is added in the group of */
/* reciprocal space rotations. */
/* Return 0 if failed */
long grg_get_reciprocal_point_group(long rec_rotations[48][3][3],
const long (*rotations)[3][3],
const long num_rot,
const long is_time_reversal)
const long is_time_reversal,
const long is_transpose)
{
long i, j, num_rot_ret, inv_exist;
const long inversion[3][3] = {
@ -266,6 +272,7 @@ long grg_get_reciprocal_point_group(long rec_rotations[48][3][3],
{ 0, 0,-1 }
};
/* Collect unique rotations */
num_rot_ret = 0;
for (i = 0; i < num_rot; i++) {
for (j = 0; j < num_rot_ret; j++) {
@ -282,8 +289,14 @@ long grg_get_reciprocal_point_group(long rec_rotations[48][3][3],
;
}
inv_exist = 0;
if (is_transpose) {
for (i = 0; i < num_rot_ret; i++) {
lagmat_transpose_matrix_l3(rec_rotations[i], rec_rotations[i]);
}
}
if (is_time_reversal) {
inv_exist = 0;
for (i = 0; i < num_rot_ret; i++) {
if (lagmat_check_identity_matrix_l3(inversion, rec_rotations[i])) {
inv_exist = 1;
@ -304,10 +317,6 @@ long grg_get_reciprocal_point_group(long rec_rotations[48][3][3],
}
}
for (i = 0; i < num_rot_ret; i++) {
lagmat_transpose_matrix_l3(rec_rotations[i], rec_rotations[i]);
}
return num_rot_ret;
err:
return 0;

3
c/grgrid.h
View File

@ -75,6 +75,7 @@ void grg_get_ir_grid_map(long *ir_grid_map,
long grg_get_reciprocal_point_group(long rec_rotations[48][3][3],
const long (*rotations)[3][3],
const long num_rot,
const long is_time_reversal);
const long is_time_reversal,
const long is_transpose);
#endif

3
c/gridsys.c
View File

@ -100,7 +100,8 @@ long gridsys_get_reciprocal_point_group(long rec_rotations[48][3][3],
return grg_get_reciprocal_point_group(rec_rotations,
rotations,
num_rot,
is_time_reversal);
is_time_reversal,
1);
}

44
c/triplet_grid.c
View File

@ -89,15 +89,16 @@ long tpk_get_ir_triplets_at_q(long *map_triplets,
const long num_rot,
const long swappable)
{
long i, num_ir;
RotMats *rec_rotations, *rotations;
long num_ir;
RotMats *rotations;
rec_rotations = bzg_alloc_RotMats(num_rot);
for (i = 0; i < num_rot; i++) {
lagmat_copy_matrix_l3(rec_rotations->mat[i], rec_rotations_in[i]);
rotations = bzg_get_reciprocal_point_group(rec_rotations_in,
num_rot,
is_time_reversal,
0);
if (rotations == NULL) {
return 0;
}
rotations = bzg_get_point_group_reciprocal(rec_rotations, is_time_reversal);
bzg_free_RotMats(rec_rotations);
num_ir = get_ir_triplets_at_q(map_triplets,
map_q,
@ -106,6 +107,7 @@ long tpk_get_ir_triplets_at_q(long *map_triplets,
rotations,
swappable);
bzg_free_RotMats(rotations);
rotations = NULL;
return num_ir;
}
@ -142,18 +144,26 @@ static long get_ir_triplets_at_q(long *map_triplets,
rot_reciprocal_q = get_point_group_reciprocal_with_q(rot_reciprocal,
D_diag,
grid_point);
num_ir_q = bzg_get_ir_grid_map(map_q,
D_diag,
PS,
rot_reciprocal_q);
grg_get_ir_grid_map(map_q,
rot_reciprocal_q->mat,
rot_reciprocal_q->size,
D_diag,
PS);
num_ir_q = 0;
for (i = 0; i < D_diag[0] * D_diag[1] * D_diag[2]; i++) {
if (map_q[i] == i) {
num_ir_q++;
}
}
if (swappable) {
num_ir_triplets = get_ir_triplets_at_q_perm_q1q2(map_triplets,
map_q,
grid_point,
D_diag,
rot_reciprocal_q,
num_ir_q);
num_ir_triplets = get_ir_triplets_at_q_perm_q1q2(map_triplets,
map_q,
grid_point,
D_diag,
rot_reciprocal_q,
num_ir_q);
} else {
num_ir_triplets = get_ir_triplets_at_q_noperm(map_triplets,
map_q,