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
« 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
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)
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)