Coverage for ase / _4 / symopt / test_fd_grad.py: 100.00%

1import numpy as np 

2import pytest 

3 

4from ase._4.symopt.relax import Relax 

5from ase.build import bulk 

6from ase.calculators.emt import EMT 

7from ase.optimize.bfgs import BFGS 

8from ase.parallel import world 

9 

10pytestmark = pytest.mark.calculator_lite 

11 

12 

13def test_generalized_coordinate_units(): 

14 """Test that displacements are in Å and stresses in eV/Å^3. 

15 

16 Even the relaxation uses generalized coordinates, they are scaled 

17 in such way, that at least in the initial configuration, moving 

18 generalized coordinate an amount x, reflects to the actual 

19 coordinate to move an amount of x. 

20 

21 |dR| = |dz(atom)| 

22 |deps| = |dz(cell)| (Cell strain) 

23 

24 """ 

25 atoms = bulk('AuAg', crystalstructure='wurtzite', a=3.24, c=5.20) 

26 print(atoms.cell.volume, 'SCALE') 

27 relax = Relax( 

28 atoms=atoms, 

29 calc=EMT, 

30 optimizer_factory=lambda atoms: BFGS(atoms, trajectory='a.traj'), 

31 symprec=0.01, 

32 comm=world, 

33 ) 

34 optimizable = relax.symmetry_adapted_atoms.__ase_optimizable__() 

35 for z in range(3): 

36 vec = np.zeros((3,)) 

37 optimizable.set_x(vec) 

38 atoms0 = optimizable.actual_atoms.copy() 

39 grad = optimizable.get_gradient() 

40 F, S = ( 

41 optimizable.actual_atoms.get_forces(), 

42 optimizable.actual_atoms.get_stress(), 

43 ) 

44 

45 print('generalized grad', grad[2], 'vs.', F) 

46 scale = np.max(np.abs(grad[2])) / np.max(np.abs(F)) 

47 print('SCALE F', scale) 

48 assert 3.9 < scale < 4.1, scale 

49 

50 print('generalized gradS S', grad[:2], 'vs', S) 

51 scale = np.max(np.abs(grad[:2])) / np.max(np.abs(S)) 

52 print('SCALE S', scale) 

53 assert 39.9 < scale < 40.1, scale 

54 

55 vec[z] = 1e-6 

56 optimizable.set_x(vec) 

57 atoms1 = optimizable.actual_atoms.copy() 

58 if z < 2: 

59 deps = np.sum((atoms1.cell - atoms0.cell) ** 2) ** 0.5 

60 print('deps', deps) 

61 else: 

62 dR = np.max( 

63 np.linalg.norm( 

64 atoms0.get_positions() - atoms1.get_positions(), axis=1 

65 ) 

66 ) 

67 print('dR', dR) 

68 

69 

70def test_fd_gradients(): 

71 atoms = bulk('AuAg', crystalstructure='wurtzite', a=3.24, c=5.20) 

72 relax = Relax( 

73 atoms=atoms, 

74 calc=EMT, 

75 optimizer_factory=lambda atoms: BFGS(atoms, trajectory='a.traj'), 

76 symprec=0.01, 

77 comm=world, 

78 ) 

79 optimizable = relax.symmetry_adapted_atoms.__ase_optimizable__() 

80 for z in range(3): 

81 vec = np.zeros((3,)) 

82 optimizable.set_x(vec) 

83 E0 = optimizable.get_value() 

84 grad = optimizable.get_gradient() 

85 vec[z] = 1e-6 

86 optimizable.set_x(vec) 

87 E1 = optimizable.get_value() 

88 print('Finite difference grad', (E1 - E0) / 1e-6) 

89 print('Gotten grad', grad[z]) 

90 div = grad[z] / ((E1 - E0) / 1e-6) 

91 assert 0.99 < div < 1.01