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
« 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
16h = 0.2
17d = 2.5
18min_vac = 4.0
19u0 = .180
20epsinf = 80.
21T = 298.15
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))
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)
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)