Coverage for gpaw/test/rpa/test_rpa_energy_N2.py: 100%

1import pytest 

2 

3from ase.units import Hartree 

4 

5from gpaw import GPAW 

6from gpaw.mpi import serial_comm 

7from gpaw.xc.rpa import RPACorrelation 

8from gpaw.hybrids.energy import non_self_consistent_energy as nsc_energy 

9 

10 

11@pytest.mark.rpa 

12@pytest.mark.response 

13def test_rpa_rpa_energy_N2(in_tmp_dir, gpw_files, scalapack): 

14 ecut = 25 

15 

16 N2_calc = GPAW(gpw_files['n2_pw'], communicator=serial_comm) 

17 

18 E_n2_pbe = N2_calc.get_potential_energy() 

19 E_n2_hf = nsc_energy(N2_calc, 'EXX') 

20 

21 rpa = RPACorrelation(N2_calc, nfrequencies=8, ecut=[ecut]) 

22 data = rpa.calculate_all_contributions() 

23 energy_w = data.energy_wi[:, 0] # we only have 1 ecut 

24 assert energy_w == pytest.approx( 

25 [-1.7415, -1.7336, -1.4659, -0.4209, -0.0295, -0.0010] 

26 + [0.0] * 2, abs=0.001) 

27 E_n2_rpa = data.energy_i * Hartree 

28 

29 N_calc = GPAW(gpw_files['n_pw'], communicator=serial_comm) 

30 E_n_pbe = N_calc.get_potential_energy() 

31 E_n_hf = nsc_energy(N_calc, 'EXX') 

32 

33 rpa = RPACorrelation(N_calc, nfrequencies=8, ecut=[ecut]) 

34 E_n_rpa = rpa.calculate() 

35 

36 print('Atomization energies:') 

37 print('PBE: ', E_n2_pbe - 2 * E_n_pbe) 

38 print('HF: ', E_n2_hf - 2 * E_n_hf) 

39 print('HF+RPA: ', E_n2_hf - 2 * E_n_hf + E_n2_rpa[0] - 2 * E_n_rpa[0]) 

40 

41 assert E_n2_rpa - 2 * E_n_rpa == pytest.approx(-1.68, abs=0.02) 

42 assert (E_n2_hf - 2 * E_n_hf) == pytest.approx( 

43 [-10.47, 3.03, 0, 0, -0.06, 2.79], abs=0.01)