Coverage for gpaw/test/response/test_mpa_poly_C.py: 100%

1import pytest 

2import numpy as np 

3from ase.units import Hartree as Ha 

4from gpaw.cgpaw import evaluate_mpa_poly as mpa_C 

5 

6 

7def mpa_py(omega, f, omegat_nGG, W_nGG, eta, factor): 

8 x1_nGG = f / (omega + omegat_nGG - 1j * eta) 

9 x2_nGG = (1.0 - f) / (omega - omegat_nGG + 1j * eta) 

10 

11 x_GG = 2 * factor * np.sum(W_nGG * (x1_nGG + x2_nGG), 

12 axis=0) 

13 

14 eps = 0.0001 / Ha 

15 xp_nGG = f / (omega + eps + omegat_nGG - 1j * eta) 

16 xp_nGG += (1.0 - f) / (omega + eps - omegat_nGG + 1j * eta) 

17 xm_nGG = f / (omega - eps + omegat_nGG - 1j * eta) 

18 xm_nGG += (1.0 - f) / (omega - eps - omegat_nGG + 1j * eta) 

19 dx_GG = 2 * factor * np.sum(W_nGG * (xp_nGG - xm_nGG) / (2 * eps), 

20 axis=0) 

21 return x_GG, dx_GG 

22 

23 

24@pytest.mark.parametrize('f', [0, 0.4, 1.0]) 

25def test_residues(f): 

26 factor = 2.0 

27 eta = 0.1 * Ha 

28 nG = 5 

29 npols = 10 

30 omega = 0.5 

31 

32 rng = np.random.default_rng(seed=1) 

33 omegat_nGG = rng.random((npols, nG, nG)) * 0.05 + 5.5 - 0.01j 

34 W_nGG = np.array(rng.random((npols, nG, nG)), dtype=complex) 

35 

36 x_GG_py, dx_GG_py = mpa_py(omega, f, omegat_nGG, W_nGG, eta, factor) 

37 

38 x_GG_C = np.empty(omegat_nGG.shape[1:], dtype=complex) 

39 dx_GG_C = np.empty(omegat_nGG.shape[1:], dtype=complex) 

40 mpa_C(x_GG_C, dx_GG_C, omega, f, omegat_nGG, W_nGG, eta, factor) 

41 

42 assert np.allclose(x_GG_py, x_GG_C, atol=1e-6)