Coverage for gpaw/test/poisson/test_pw_charged.py: 100%

30 statements  

« prev     ^ index     » next       coverage.py v7.7.1, created at 2025-07-14 00:18 +0000

1import pytest 

2from ase import Atoms 

3from gpaw import GPAW 

4import numpy as np 

5from gpaw.grid_descriptor import GridDescriptor 

6from gpaw.pw.descriptor import PWDescriptor 

7from gpaw.pw.poisson import ChargedReciprocalSpacePoissonSolver as CRSPC 

8 

9 

10def test_charged_pw_poisson(): 

11 """Calculate Coulomb energy of Gaussian charge and compare to analytic 

12 result. 

13 """ 

14 n = 50 

15 L = n / 40 * 8.0 

16 charge = 1.5 

17 gd = GridDescriptor((n, n, n), (L, L, L)) 

18 pd = PWDescriptor(None, gd) 

19 ps = CRSPC(pd, charge, alpha=1.1) 

20 a = 0.8 

21 C = gd.cell_cv.sum(0) / 2 

22 rho = -np.exp(-1 / (4 * a) * pd.G2_qG[0] + 

23 1j * (pd.get_reciprocal_vectors() @ C)) / gd.dv 

24 v = np.empty_like(rho) 

25 e = ps._solve(v, rho * charge) 

26 w = gd.collect(pd.ifft(v), broadcast=True) 

27 for pot, d in [(w[n // 2, n // 2, 0], L / 2), 

28 (w[n // 2, 0, 0], L / 2 * 2**0.5), 

29 (w[0, 0, 0], L / 2 * 3**0.5)]: 

30 assert pot == pytest.approx(-charge / d, abs=0.002) 

31 assert e == pytest.approx((a / 2 / np.pi)**0.5 * charge**2) 

32 

33 

34def test_pw_proton(): 

35 """Check that the energy of a proton is 0.0.""" 

36 proton = Atoms('H') 

37 proton.center(vacuum=2.0) 

38 proton.calc = GPAW(mode='pw', charge=1) 

39 e = proton.get_potential_energy() 

40 e += proton.calc.get_reference_energy() 

41 assert e == pytest.approx(0.0, abs=0.004)