Coverage for gpaw/test/xc/test_pplda.py: 100%

1from ase import Atoms 

2from gpaw import GPAW 

3import numpy as np 

4 

5# spin paired H2 

6 

7 

8def test_xc_pplda(): 

9 d = 0.75 

10 h2 = Atoms('H2', [[0, 0, 0], [0, 0, d]]) 

11 h2.center(vacuum=2.) 

12 

13 e = np.array([]) 

14 f = np.array([]) 

15 

16 for xc in ['LDA', 'PPLDA']: 

17 calc = GPAW(nbands=-1, 

18 xc=xc, 

19 convergence={'eigenstates': 1.e-9}, 

20 txt=None, 

21 mode='lcao', 

22 basis='szp(dzp)', 

23 h=0.3) 

24 h2.calc = calc 

25 e = np.append(e, h2.get_potential_energy()) 

26 f = np.append(f, h2.get_forces()) 

27 del calc 

28 

29 print(e[0], e[1], np.abs(e[0] - e[1])) 

30 print(f[0], f[1], np.sum(np.abs(f[0] - f[1]))) 

31 assert np.abs(e[0] - e[1]) < 1.e-8 

32 assert np.sum(np.abs(f[0] - f[1])) < 1.e-8 

33 

34 # spin polarized O2 

35 d = 1.2 

36 o2 = Atoms('O2', [[0, 0, 0], [0, 0, d]], magmoms=[1., 1.]) 

37 o2.center(vacuum=2.) 

38 

39 e = np.array([]) 

40 f = np.array([]) 

41 

42 for xc in ['LDA', 'PPLDA']: 

43 calc = GPAW(nbands=-2, xc=xc, convergence={'energy': 1.e-9}, txt=None, 

44 mode='lcao', basis='szp(dzp)', h=0.3) 

45 o2.calc = calc 

46 e = np.append(e, o2.get_potential_energy()) 

47 f = np.append(f, o2.get_forces()) 

48 del calc 

49 

50 print(e[0], e[1], np.abs(e[0] - e[1])) 

51 print(f[0], f[1], np.sum(np.abs(f[0] - f[1]))) 

52 assert np.abs(e[0] - e[1]) < 1.e-8 

53 assert np.sum(np.abs(f[0] - f[1])) < 1.e-4