Add C gridsys test_gridsys_get_triplets_at_q_wurtzite_force_SNF

This commit is contained in:
Atsushi Togo 2023-01-26 15:39:48 +09:00
parent b326e5867a
commit 3f29c73c7f
5 changed files with 528 additions and 5 deletions

View File

@ -26,11 +26,13 @@ if(NOT GTest_FOUND)
endif()
foreach(testcase IN ITEMS test_gridsys test_niggli)
add_executable(${testcase} ${testcase}.cpp)
add_executable(${testcase} ${CMAKE_CURRENT_SOURCE_DIR}/${testcase}.cpp
${CMAKE_CURRENT_SOURCE_DIR}/utils.c)
target_link_libraries(
${testcase}
PUBLIC gridsys
PRIVATE GTest::gtest GTest::gtest_main)
target_include_directories(${testcase} PUBLIC ${PROJECT_SOURCE_DIR}/c)
target_include_directories(${testcase} PUBLIC ${PROJECT_SOURCE_DIR}/c
${CMAKE_CURRENT_SOURCE_DIR})
gtest_discover_tests(${testcase})
endforeach()

View File

@ -61,6 +61,22 @@ const long wurtzite_rec_rotations_without_time_reversal[12][3][3] = {
{{1, 0, 0}, {-1, -1, 0}, {0, 0, 1}}, {{0, -1, 0}, {-1, 0, 0}, {0, 0, 1}},
{{-1, -1, 0}, {0, 1, 0}, {0, 0, 1}}, {{-1, 0, 0}, {1, 1, 0}, {0, 0, 1}}};
// Transformed point group operations of wurtzite {R^T} (without time reversal)
// D_diag=[1, 5, 15], Q=[[-1, 0, -6], [0, -1, 0], [-1, 0, -5]]
const long wurtzite_tilde_rec_rotations_without_time_reversal[12][3][3] = {
{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
{{1, -1, 0}, {-5, 0, -2}, {0, 3, 1}},
{{6, -1, 2}, {-5, -1, -2}, {-15, 3, -5}},
{{11, 0, 4}, {0, -1, 0}, {-30, 0, -11}},
{{11, 1, 4}, {5, 0, 2}, {-30, -3, -11}},
{{6, 1, 2}, {5, 1, 2}, {-15, -3, -5}},
{{6, -1, 2}, {5, 0, 2}, {-15, 3, -5}},
{{1, -1, 0}, {0, -1, 0}, {0, 3, 1}},
{{1, 0, 0}, {-5, -1, -2}, {0, 0, 1}},
{{6, 1, 2}, {-5, 0, -2}, {-15, -3, -5}},
{{11, 1, 4}, {0, 1, 0}, {-30, -3, -11}},
{{11, 0, 4}, {5, 1, 2}, {-30, 0, -11}}};
// Symmetry operations of wurtzite 1x1x2 {R}
const long wurtzite112_symmetry_operations[24][3][3] = {
{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}, {{1, -1, 0}, {1, 0, 0}, {0, 0, 1}},
@ -92,6 +108,20 @@ const long AgNO2_tilde_rec_rotations[8][3][3] = {
{{-1, 0, 0}, {0, -11, -2}, {0, 60, 11}},
{{-5, 12, 2}, {-12, 35, 6}, {60, -180, -31}}};
const long AgNO2_tilde_rec_rotations_with_time_reversal_mesh12[8][3][3] = {
{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}, {{1, 0, 0}, {0, -1, 0}, {4, 0, -1}},
{{1, 0, 0}, {0, 1, 0}, {4, 4, -1}}, {{1, 0, 0}, {0, -1, 0}, {0, -4, 1}},
{{-1, 0, 0}, {0, -1, 0}, {0, 0, -1}}, {{-1, 0, 0}, {0, 1, 0}, {-4, 0, 1}},
{{-1, 0, 0}, {0, -1, 0}, {-4, -4, 1}}, {{-1, 0, 0}, {0, 1, 0}, {0, 4, -1}},
};
// D_diag=[2, 2, 8], Q=[[0, 0, 1], [1, 0, -1], [0, 1, -1]]
const long AgNO2_tilde_rec_rotations_without_time_reversal_mesh12[4][3][3] = {
{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
{{1, 0, 0}, {0, -1, 0}, {4, 0, -1}},
{{1, 0, 0}, {0, 1, 0}, {4, 4, -1}},
{{1, 0, 0}, {0, -1, 0}, {0, -4, 1}}};
/**
* @brief gridsys_get_all_grid_addresses test
* Return all GR-grid addresses of {(x, y, z)} where x runs fastest.
@ -898,3 +928,100 @@ TEST(test_gridsys, test_gridsys_get_triplets_at_q_wurtzite) {
}
}
}
/**
* @brief gridsys_get_triplets_at_q by AgNO2 rotations
* @details Four patterns, is_time_reversal x swappable, are tested.
*/
TEST(test_gridsys, test_gridsys_get_triplets_at_q_AgNO2) {
long D_diag[3] = {2, 2, 8};
long grid_point = 1;
long map_triplets[32], map_q[32];
long i, j, k, num_triplets;
long ref_num_triplets[2][2] = {{8, 16}, {12, 24}};
long is_time_reversal[2] = {1, 0};
long swappable[2] = {1, 0};
long count = 0;
long ref_map_triplets[4][32] = {
{0, 0, 2, 2, 4, 4, 6, 6, 8, 8, 10, 10, 12, 12, 6, 6,
16, 16, 2, 2, 12, 12, 6, 6, 8, 8, 10, 10, 4, 4, 6, 6},
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 6, 7,
16, 17, 2, 3, 12, 13, 6, 7, 8, 9, 10, 11, 4, 5, 6, 7},
{0, 0, 2, 2, 4, 5, 6, 7, 8, 9, 10, 10, 12, 13, 7, 6,
16, 16, 2, 2, 13, 12, 6, 7, 9, 8, 10, 10, 5, 4, 7, 6},
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 2, 3, 20, 21, 6, 7, 24, 25, 10, 11, 28, 29, 14, 15}};
long ref_map_q[4][32] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 6, 7,
16, 17, 2, 3, 12, 13, 6, 7, 8, 9, 10, 11, 4, 5, 6, 7},
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 6, 7,
16, 17, 2, 3, 12, 13, 6, 7, 8, 9, 10, 11, 4, 5, 6, 7},
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 2, 3, 20, 21, 6, 7, 24, 25, 10, 11, 28, 29, 14, 15},
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 2, 3, 20, 21, 6, 7, 24, 25, 10, 11, 28, 29, 14, 15}};
for (i = 0; i < 2; i++) {
for (j = 0; j < 2; j++) {
num_triplets = gridsys_get_triplets_at_q(
map_triplets, map_q, grid_point, D_diag, is_time_reversal[i], 4,
AgNO2_tilde_rec_rotations_without_time_reversal_mesh12,
swappable[j]);
ASSERT_EQ(ref_num_triplets[i][j], num_triplets);
for (k = 0; k < 32; k++) {
ASSERT_EQ(ref_map_triplets[count][k], map_triplets[k]);
ASSERT_EQ(ref_map_q[count][k], map_q[k]);
}
count++;
}
}
}
/**
* @brief gridsys_get_triplets_at_q by wurtzite rotations with and without
* force_SNF (i.e., transformed or not transformed rotations)
* @details Four patterns, is_time_reversal x swappable, are tested.
*/
TEST(test_gridsys, test_gridsys_get_triplets_at_q_wurtzite_force_SNF) {
long D_diag[2][3] = {{1, 5, 15}, {5, 5, 3}};
long grid_point = 1;
long map_triplets[75], map_q[75];
long i, j, k, ll, num_triplets;
long ref_unique_elems[4][2] = {{18, 30}, {24, 45}, {30, 30}, {45, 45}};
long is_time_reversal[2] = {1, 0};
long swappable[2] = {1, 0};
long rec_rotations[2][12][3][3];
for (i = 0; i < 2; i++) {
for (j = 0; j < 12; j++) {
for (k = 0; k < 3; k++) {
for (ll = 0; ll < 3; ll++) {
if (i == 0) {
rec_rotations[i][j][k][ll] =
wurtzite_tilde_rec_rotations_without_time_reversal
[j][k][ll];
} else {
rec_rotations[i][j][k][ll] =
wurtzite_rec_rotations_without_time_reversal[j][k]
[ll];
}
}
}
}
}
for (i = 0; i < 2; i++) { // force_SNF True or False
for (j = 0; j < 2; j++) { // swappable
for (k = 0; k < 2; k++) { // is_time_reversal
num_triplets = gridsys_get_triplets_at_q(
map_triplets, map_q, grid_point, D_diag[i],
is_time_reversal[k], 12, rec_rotations[i], swappable[j]);
ASSERT_EQ(ref_unique_elems[j * 2 + k][0], num_triplets);
ASSERT_EQ(ref_unique_elems[j * 2 + k][0],
get_num_unique_elems(map_triplets, 75));
ASSERT_EQ(ref_unique_elems[j * 2 + k][1],
get_num_unique_elems(map_q, 75));
}
}
}
}

View File

@ -1,4 +1,7 @@
#include "utils.h"
#include <stdio.h>
#include <stdlib.h>
void lagmat_multiply_matrix_vector_l3(long v[3], const long a[3][3],
const long b[3]) {
@ -45,3 +48,25 @@ long lagmat_check_identity_matrix_l3(const long a[3][3], const long b[3][3]) {
return 1;
}
}
long get_num_unique_elems(const long array[], const long array_size) {
long i, num_unique_elems;
long *unique_elems;
num_unique_elems = 0;
unique_elems = (long *)malloc(sizeof(long) * array_size);
for (i = 0; i < array_size; i++) {
unique_elems[i] = 0;
}
for (i = 0; i < array_size; i++) {
unique_elems[array[i]]++;
}
for (i = 0; i < array_size; i++) {
if (unique_elems[i]) {
num_unique_elems++;
}
}
free(unique_elems);
unique_elems = NULL;
return num_unique_elems;
}

View File

@ -7,4 +7,5 @@ void lagmat_multiply_matrix_l3(long m[3][3], const long a[3][3],
const long b[3][3]);
void lagmat_copy_matrix_l3(long a[3][3], const long b[3][3]);
long lagmat_check_identity_matrix_l3(const long a[3][3], const long b[3][3]);
long get_num_unique_elems(const long array[], const long array_size);
#endif // __test_utils_H__

View File

@ -1,6 +1,7 @@
"""Test for triplets.py."""
import numpy as np
import pytest
from phonopy import Phonopy
from phonopy.structure.atoms import PhonopyAtoms
from phonopy.structure.symmetry import Symmetry
@ -492,8 +493,375 @@ def test_get_triplets_reciprocal_mesh_at_q(aln_cell: PhonopyAtoms, params):
is_time_reversal=params[1],
swappable=params[0],
)
print(",".join(["%d" % x for x in map_triplets]))
print(",".join(["%d" % x for x in map_q]))
print(len(np.unique(map_triplets)))
# print(",".join(["%d" % x for x in map_triplets]))
# print(",".join(["%d" % x for x in map_q]))
# print(len(np.unique(map_triplets)))
np.testing.assert_equal(ref_map_triplets[params[2]], map_triplets)
np.testing.assert_equal(ref_map_q[params[2]], map_q)
@pytest.mark.parametrize(
"params",
[(True, True, 0), (False, True, 1), (True, False, 2), (False, False, 3)],
)
def test_get_triplets_reciprocal_mesh_at_q_agno2(agno2_cell: PhonopyAtoms, params):
"""Test BZGrid with shift using AgNO2."""
ref_map_triplets = [
[
0,
0,
2,
2,
4,
4,
6,
6,
8,
8,
10,
10,
12,
12,
6,
6,
16,
16,
2,
2,
12,
12,
6,
6,
8,
8,
10,
10,
4,
4,
6,
6,
],
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
6,
7,
16,
17,
2,
3,
12,
13,
6,
7,
8,
9,
10,
11,
4,
5,
6,
7,
],
[
0,
0,
2,
2,
4,
5,
6,
7,
8,
9,
10,
10,
12,
13,
7,
6,
16,
16,
2,
2,
13,
12,
6,
7,
9,
8,
10,
10,
5,
4,
7,
6,
],
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
14,
15,
16,
17,
2,
3,
20,
21,
6,
7,
24,
25,
10,
11,
28,
29,
14,
15,
],
]
ref_map_q = [
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
6,
7,
16,
17,
2,
3,
12,
13,
6,
7,
8,
9,
10,
11,
4,
5,
6,
7,
],
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
6,
7,
16,
17,
2,
3,
12,
13,
6,
7,
8,
9,
10,
11,
4,
5,
6,
7,
],
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
14,
15,
16,
17,
2,
3,
20,
21,
6,
7,
24,
25,
10,
11,
28,
29,
14,
15,
],
[
0,
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
14,
15,
16,
17,
2,
3,
20,
21,
6,
7,
24,
25,
10,
11,
28,
29,
14,
15,
],
]
grid_point = 1
mesh = 12
ph = Phonopy(agno2_cell, supercell_matrix=[1, 1, 1], primitive_matrix="auto")
bzgrid = BZGrid(
mesh,
lattice=ph.primitive.cell,
symmetry_dataset=ph.primitive_symmetry.dataset,
use_grg=True,
is_time_reversal=False,
)
np.testing.assert_equal([2, 2, 8], bzgrid.D_diag)
np.testing.assert_equal([[0, 0, 1], [1, 0, -1], [0, 1, -1]], bzgrid.Q)
map_triplets, map_q = _get_triplets_reciprocal_mesh_at_q(
grid_point,
bzgrid.D_diag,
bzgrid.rotations,
is_time_reversal=params[1],
swappable=params[0],
)
# print(",".join(["%d" % x for x in map_triplets]))
# print(",".join(["%d" % x for x in map_q]))
np.testing.assert_equal(ref_map_triplets[params[2]], map_triplets)
np.testing.assert_equal(ref_map_q[params[2]], map_q)
@pytest.mark.parametrize(
"params",
[ # force_SNF, swappable, is_time_reversal
(True, True, True, 0),
(True, True, False, 1),
(True, False, True, 2),
(True, False, False, 3),
(False, True, True, 4),
(False, True, False, 5),
(False, False, True, 6),
(False, False, False, 7),
],
)
def test_get_triplets_reciprocal_mesh_at_q_wurtzite_force(
aln_cell: PhonopyAtoms, params
):
"""Test BZGrid with shift using wurtzite with and without force_SNF."""
grid_point = 1
mesh = 14
ph = Phonopy(aln_cell, supercell_matrix=[1, 1, 1], primitive_matrix="auto")
bzgrid = BZGrid(
mesh,
lattice=ph.primitive.cell,
symmetry_dataset=ph.primitive_symmetry.dataset,
use_grg=True,
force_SNF=params[0],
is_time_reversal=False,
)
for r in bzgrid.rotations:
print("{")
for v in r:
print("{%d, %d, %d}," % tuple(v))
print("},")
ref_unique_elems = [[18, 30], [24, 45], [30, 30], [45, 45]]
if params[0]:
np.testing.assert_equal([1, 5, 15], bzgrid.D_diag)
np.testing.assert_equal([[-1, 0, -6], [0, -1, 0], [-1, 0, -5]], bzgrid.Q)
else:
np.testing.assert_equal([5, 5, 3], bzgrid.D_diag)
np.testing.assert_equal(np.eye(3, dtype=int), bzgrid.Q)
map_triplets, map_q = _get_triplets_reciprocal_mesh_at_q(
grid_point,
bzgrid.D_diag,
bzgrid.rotations,
is_time_reversal=params[2],
swappable=params[1],
)
# "% 4" means that expectation of the same values with and without force_SNF.
np.testing.assert_equal(
ref_unique_elems[params[3] % 4],
[len(np.unique(map_triplets)), len(np.unique(map_q))],
)