Coverage for gpaw/test/fd_ops/test_laplace.py: 100%

1import pytest 

2from gpaw.mpi import world 

3from math import sin, cos, pi 

4import numpy as np 

5from gpaw.fd_operators import GUCLaplace as Laplace 

6from gpaw.fd_operators import Gradient 

7from gpaw.grid_descriptor import GridDescriptor 

8 

9pytestmark = pytest.mark.skipif(world.size > 1, 

10 reason='world.size > 1') 

11 

12 

13def test_fd_ops_laplace(): 

14 cells = [ 

15 ('distorted hexagonal', 4, 

16 [(1, 0, 0), 

17 (1.02 * cos(pi / 3 - 0.02), 1.02 * sin(pi / 3 - 0.02), 0), 

18 (0, 0, 1.0)]), 

19 ('hexagonal', 4, 

20 [(1, 0, 0), 

21 (0.5, 3**0.5 / 2, 0), 

22 (0, 0, 1.1)]), 

23 ('fcc', 6, 

24 [(0, 1, 1), 

25 (1, 0, 1), 

26 (1, 1, 0)]), 

27 ('fcc-alternative', 6, 

28 [(1, 0, 0), 

29 (0.5, 3**0.5 / 2, 0), 

30 (0.5, 3**0.5 / 6, (2 / 3)**0.5)]), 

31 ('bcc', 4, 

32 [(-1, 1, 1), 

33 (1, -1, 1), 

34 (1, 1, -1)]), 

35 ('sc', 3, 

36 [1.1, 1.02, 1.03]), 

37 ('distorted sc', 6, 

38 [(1, 0, 0), 

39 (0.01, 1, 0), 

40 (0, 0.02, 1)]), 

41 ('rocksalt', 6, 

42 [(2 * np.sqrt(1 / 3), np.sqrt(1 / 8), -np.sqrt(1 / 24)), 

43 (2 * np.sqrt(1 / 3), -np.sqrt(1 / 8), -np.sqrt(1 / 24)), 

44 (2 * np.sqrt(1 / 3), 0, np.sqrt(1 / 6))]), 

45 ('nasty', 6, 

46 [(1, 0, 0), 

47 (0.0001, 1.03, 0), 

48 (0.0001, 0.0001, 1.0)]), 

49 ('Mike', 6, 

50 5 * np.array([(5.565 / 28, 0, 0), 

51 (0.0001 / 28, 5.565 / 28, 0), 

52 (0.0001 / 24, 0.0001 / 24, 4.684 / 24)])), 

53 ('MnO', 6, 

54 [(1, 0.5, 0.5), (0.5, 1, 0.5), (0.5, 0.5, 1)]), 

55 ('Jelver', 6, 

56 [[6.658433626216136, 0.1711724130951983, -0.04300038284455887], 

57 [3.4774564712755938, 5.843379292501022, 0.01599293966594096], 

58 [-0.10777038306906983 * 0.43, 0.10850460815311265 * 0.43, 

59 15.26098014321118 * 0.43]])] 

60 

61 for name, D, cell in cells: 

62 if name == 'Jelver': 

63 # Strange one! 

64 continue 

65 

66 print('------------------') 

67 print(name, D) 

68 print(cell[0]) 

69 print(cell[1]) 

70 print(cell[2]) 

71 for n in range(1, 5): 

72 N = 2 * n + 2 

73 gd = GridDescriptor((N, N, N), cell) 

74 b_g = gd.zeros() 

75 r_gv = gd.get_grid_point_coordinates().transpose((1, 2, 3, 0)) 

76 c_v = gd.cell_cv.sum(0) / 2 

77 r_gv -= c_v 

78 lap = Laplace(gd, n=n) 

79 grad_v = [Gradient(gd, v, n=n) for v in range(3)] 

80 assert lap.npoints == D * 2 * n + 1 

81 print([g.npoints for g in grad_v], n, D) 

82 for m in range(0, 2 * n + 1): 

83 for ix in range(m + 1): 

84 for iy in range(m - ix + 1): 

85 iz = m - ix - iy 

86 a_g = (r_gv**(ix, iy, iz)).prod(3) 

87 if ix + iy + iz == 2 and max(ix, iy, iz) == 2: 

88 r = 2.0 

89 else: 

90 r = 0.0 

91 lap.apply(a_g, b_g) 

92 e = b_g[n + 1, n + 1, n + 1] - r 

93 assert abs(e) < 2e-12, e 

94 for v in range(3): 

95 grad_v[v].apply(a_g, b_g) 

96 if m == 1 and [ix, iy, iz][v] == 1: 

97 r = 1 

98 else: 

99 r = 0 

100 e = b_g[n + 1, n + 1, n + 1] - r 

101 assert abs(e) < 4e-12, (n, ix, iy, iz, r, v, e)