1import pytest
2import numpy as np
3from phonopy.harmonic.force_constants import cutoff_force_constants
4from phonopy.structure.cells import get_primitive
5
6
7fc_1_10_ref = [-0.037549, 0.000000, 0.000000,
8               0.000000, 0.002415, -0.001746,
9               0.000000, -0.001746, 0.002415]
10
11fc_1_10_nofcsym_ref = [-0.005051, 0.000000, 0.000000,
12                       0.000000, 0.094457, 0.000000,
13                       0.000000, 0.000000, -0.020424]
14
15fc_1_10_compact_fcsym_ref = [-0.004481, 0.000000, 0.000000,
16                             0.000000, 0.095230, 0.000000,
17                             0.000000, 0.000000, -0.019893]
18
19
20def test_fc(ph_nacl):
21    fc_1_10 = ph_nacl.force_constants[1, 10]
22    # print("".join(["%f, " % v for v in fc_1_10.ravel()]))
23    np.testing.assert_allclose(fc_1_10.ravel(), fc_1_10_ref, atol=1e-5)
24
25
26def test_fc_nofcsym(ph_nacl_nofcsym):
27    fc_1_10 = ph_nacl_nofcsym.force_constants[1, 10]
28    # print("".join(["%f, " % v for v in fc_1_10.ravel()]))
29    np.testing.assert_allclose(fc_1_10.ravel(), fc_1_10_nofcsym_ref, atol=1e-5)
30
31
32def test_fc_compact_fcsym(ph_nacl_compact_fcsym):
33    fc_1_10 = ph_nacl_compact_fcsym.force_constants[1, 10]
34    # print("".join(["%f, " % v for v in fc_1_10.ravel()]))
35    np.testing.assert_allclose(
36        fc_1_10.ravel(), fc_1_10_compact_fcsym_ref, atol=1e-5)
37
38
39@pytest.mark.parametrize("is_compact", [False, True])
40def test_fc_cutoff_radius(ph_nacl, ph_nacl_compact_fcsym, is_compact):
41    if is_compact:
42        ph = ph_nacl_compact_fcsym
43    else:
44        ph = ph_nacl
45
46    # Need restore fc because fc are overwritten.
47    fc_orig = ph.force_constants.copy()
48    ph.set_force_constants_zero_with_radius(4.0)
49    changed = (np.abs(fc_orig - ph.force_constants) > 1e-8)
50    ph.force_constants = fc_orig
51
52    if is_compact:
53        assert np.sum(changed) == 534
54    else:
55        assert np.sum(changed) == 17088
56
57
58@pytest.mark.parametrize("is_compact,store_dense_svecs",
59                         [(False, False), (False, True),
60                          (True, False), (True, True)])
61def test_fc_cutoff_radius_svecs(ph_nacl, ph_nacl_compact_fcsym,
62                                is_compact, store_dense_svecs):
63    if is_compact:
64        ph = ph_nacl_compact_fcsym
65    else:
66        ph = ph_nacl
67
68    fc = ph.force_constants.copy()
69    primitive_matrix = np.dot(np.linalg.inv(ph.supercell_matrix),
70                              ph.primitive_matrix)
71    primitive = get_primitive(ph.supercell, primitive_matrix,
72                              store_dense_svecs=store_dense_svecs)
73
74    cutoff_force_constants(fc, ph.supercell, primitive, 4.0)
75    changed = (np.abs(ph.force_constants - fc) > 1e-8)
76
77    if is_compact:
78        assert np.sum(changed) == 534
79    else:
80        assert np.sum(changed) == 17088
81