Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

""" 

This module contains various support/utility functions. 

""" 

 

import numpy as np 

from ase.geometry import find_mic 

from .io.logging import logger 

 

 

logger = logger.getChild('utilities') 

 

 

def get_displacements(atoms, atoms_ideal): 

""" 

Computes the smallest possible displacements from positions and 

ideal positions given a cell metric. 

 

Notes 

----- 

* uses :func:`ase.geometry.find_mic` 

* assumes `pbc=[True, True, True]`. 

 

Parameters 

---------- 

atoms : ase.Atoms 

configuration with displaced atoms 

atoms_ideal : ase.Atoms 

ideal configuration relative to which displacements are computed 

 

Returns 

------- 

numpy.ndarray 

wrapped displacements 

""" 

35 ↛ 36line 35 didn't jump to line 36, because the condition on line 35 was never true if not np.array_equal(atoms.numbers, atoms_ideal.numbers): 

raise ValueError('Atomic numbers do not match.') 

 

displacements = [] 

for pos, ideal_pos in zip(atoms.positions, atoms_ideal.positions): 

v_ij = np.array([pos - ideal_pos]) 

displacements.append(find_mic(v_ij, atoms_ideal.cell, pbc=True)[0][0]) 

return np.array(displacements) 

 

 

def prepare_structures(structures, atoms_ideal, calc): 

""" 

Prepares a set of structures in the format suitable for a 

:class:`StructureContainer <hiphive.StructureContainer>`. 

 

**Note:** Changes the structures in place. 

""" 

for atoms in structures: 

atoms.set_calculator(calc) 

forces = atoms.get_forces() 

displacements = get_displacements(atoms, atoms_ideal) 

atoms.positions = atoms_ideal.get_positions() 

atoms.new_array('forces', forces) 

atoms.new_array('displacements', displacements) 

 

 

def get_neighbor_shells(atoms, cutoff, dist_tol=1e-5): 

""" Gets list of neighbor shells. 

 

Distances are grouped into shells via the following algorithm: 

 

1. Find smallest atomic distance `d_min` 

 

2. Find all pair distances in the range `d_min + 1 * dist_tol` 

 

3. Construct a shell from these and pop them from distance list 

 

4. Go to 1. 

 

Parameters 

---------- 

atoms : ase.Atoms 

Atoms used for finding shells 

cutoff : float 

exclude neighbor shells which have a distance larger than cutoff 

dist_tol : float 

distance tolerance 

 

Returns 

------- 

list(Shell) 

neighbor shells 

""" 

 

# TODO: Remove this once this feature have been in ASE long enough 

try: 

from ase.neighborlist import neighbor_list 

except ImportError: 

raise ImportError('get_neighbor_shells uses new functionality from ASE' 

', please update ASE in order to use this function.') 

 

# get distances 

ijd = neighbor_list('ijd', atoms, cutoff) 

ijd = list(zip(*ijd)) 

ijd.sort(key=lambda x: x[2]) 

 

# sort into shells 

symbols = atoms.get_chemical_symbols() 

shells = [] 

for i, j, d in ijd: 

types = tuple(sorted([symbols[i], symbols[j]])) 

for shell in shells: 

if abs(d - shell.distance) < dist_tol and types == shell.types: 

shell.count += 1 

break 

else: 

shell = Shell(types, d, 1) 

shells.append(shell) 

shells.sort(key=lambda x: (x.distance, x.types, x.count)) 

 

# warning if two shells are within 2 * tol 

for i, s1 in enumerate(shells): 

for j, s2 in enumerate(shells[i+1:]): 

if s1.types != s2.types: 

continue 

if not s1.distance < s2.distance - 2 * dist_tol: 

logger.warning('Found two shells within 2 * dist_tol') 

 

return shells 

 

 

class Shell: 

""" 

Neighbor Shell class 

 

Parameters 

---------- 

types : list or tuple 

atomic types for neighbor shell 

distance : float 

interatomic distance for neighbor shell 

count : int 

number of pairs in the neighbor shell 

""" 

 

def __init__(self, types, distance, count=0): 

self.types = types 

self.distance = distance 

self.count = count 

 

def __str__(self): 

s = '{}-{} distance: {:10.6f} count: {}'.format( 

*self.types, self.distance, self.count) 

return s 

 

__repr__ = __str__