Coverage for gpaw/test/solvation/test_forces_symmetry.py: 100%

24 statements  

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

1from ase import Atoms 

2import pytest 

3from ase.units import Pascal, m 

4from gpaw.solvation import ( 

5 SolvationGPAW, 

6 EffectivePotentialCavity, 

7 Power12Potential, 

8 LinearDielectric, 

9 KB51Volume, 

10 GradientSurface, 

11 VolumeInteraction, 

12 SurfaceInteraction, 

13 LeakedDensityInteraction) 

14import numpy as np 

15 

16h = 0.2 

17d = 2.5 

18min_vac = 4.0 

19u0 = .180 

20epsinf = 80. 

21T = 298.15 

22 

23 

24def test_solvation_forces_symmetry(): 

25 xy_cell = np.ceil((min_vac * 2.) / h / 8.) * 8. * h 

26 z_cell = np.ceil((min_vac * 2. + d) / h / 8.) * 8. * h 

27 atoms = Atoms( 

28 'NaCl', positions=( 

29 (xy_cell / 2., xy_cell / 2., z_cell / 2. - d / 2.), 

30 (xy_cell / 2., xy_cell / 2., z_cell / 2. + d / 2.) 

31 ) 

32 ) 

33 atoms.set_cell((xy_cell, xy_cell, z_cell)) 

34 

35 atoms.calc = SolvationGPAW( 

36 mode='fd', 

37 xc='PBE', 

38 h=h, 

39 setups={'Na': '1'}, 

40 cavity=EffectivePotentialCavity( 

41 effective_potential=Power12Potential(u0=u0), 

42 temperature=T, 

43 volume_calculator=KB51Volume(), 

44 surface_calculator=GradientSurface()), 

45 dielectric=LinearDielectric(epsinf=epsinf), 

46 # parameters chosen to give ~ 1eV for each interaction 

47 interactions=[ 

48 VolumeInteraction(pressure=-1e9 * Pascal), 

49 SurfaceInteraction(surface_tension=100. * 1e-3 * Pascal * m), 

50 LeakedDensityInteraction(voltage=10.)]) 

51 F = atoms.calc.get_forces(atoms) 

52 

53 difference = F[0][2] + F[1][2] 

54 print(difference) 

55 assert difference == pytest.approx(.0, abs=0.02) 

56 F[0][2] = F[1][2] = 0.0 

57 print(np.abs(F)) 

58 assert np.abs(F) == pytest.approx(.0, abs=1e-10)