 r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1"""Implements the dimensionality scoring parameter.

3Method is described in:

4Definition of a scoring parameter to identify low-dimensional materials

5components

6P.M. Larsen, M. Pandey, M. Strange, and K. W. Jacobsen

7Phys. Rev. Materials 3 034003, 2019

8https://doi.org/10.1103/PhysRevMaterials.3.034003

9"""

10import numpy as np

11from collections import namedtuple

12from ase.geometry.dimensionality import rank_determination

13from ase.geometry.dimensionality import topology_scaling

14from ase.geometry.dimensionality.bond_generator import next_bond

17KInterval = namedtuple('KInterval', 'dimtype score a b h components cdim')

20def f(x):

21 if x == float("inf"):

22 return 1

23 k = 1 / 0.15**2

24 return k * max(0, x - 1)**2 / (1. + k * max(0, x - 1)**2)

27def calculate_score(a, b):

28 return f(b) - f(a)

31def reduced_histogram(h):

32 h = [int(e > 0) for e in h]

33 return tuple(h)

36def build_dimtype(h):

37 h = reduced_histogram(h)

38 return ''.join([str(i) for i, e in enumerate(h) if e > 0]) + 'D'

41def build_kinterval(a, b, h, components, cdim, score=None):

42 if score is None:

43 score = calculate_score(a, b)

45 return KInterval(dimtype=build_dimtype(h), score=score,

46 a=a, b=b, h=h, components=components, cdim=cdim)

49def merge_intervals(intervals):

50 """Merges intervals of the same dimensionality type.

52 For example, two histograms with component histograms [10, 4, 0, 0] and

53 [6, 2, 0, 0] are both 01D structures so they will be merged.

55 Intervals are merged by summing the scores, and taking the minimum and

56 maximum k-values. Component IDs in the merged interval are taken from the

57 interval with the highest score.

59 On rare occasions, intervals to be merged are not adjacent. In this case,

60 the score of the merged interval is not equal to the score which would be

61 calculated from its k-interval. This is necessary to maintain the property

62 that the scores sum to 1.

63 """

64 dimtypes = set([e.dimtype for e in intervals])

66 merged_intervals = []

67 for dimtype in dimtypes:

68 relevant = [e for e in intervals if e.dimtype == dimtype]

69 combined_score = sum([e.score for e in relevant])

70 amin = min([e.a for e in relevant])

71 bmax = max([e.b for e in relevant])

72 best = max(relevant, key=lambda x: x.score)

73 merged = build_kinterval(amin, bmax, best.h, best.components,

74 best.cdim, score=combined_score)

75 merged_intervals.append(merged)

76 return merged_intervals

79def build_kintervals(atoms, method_name):

80 """The interval analysis is performed by inserting bonds one at a time

81 until the component analysis finds a single component."""

82 method = {'RDA': rank_determination.RDA,

83 'TSA': topology_scaling.TSA}[method_name]

85 assert all([e in [0, 1] for e in atoms.pbc])

86 num_atoms = len(atoms)

87 intervals = []

88 kprev = 0

89 calc = method(num_atoms)

90 hprev = calc.check()

91 components_prev, cdim_prev = calc.get_components()

93 """

94 The end state is a single component, whose dimensionality depends on

95 the periodic boundary conditions:

96 """

97 end_state = np.zeros(4)

98 end_dim = sum(atoms.pbc)

99 end_state[end_dim] = 1

100 end_state = tuple(end_state)

102 # Insert each new bond into the component graph.

103 for (k, i, j, offset) in next_bond(atoms):

104 calc.insert_bond(i, j, offset)

105 h = calc.check()

106 if h == hprev: # Test if any components were merged

107 continue

109 components, cdim = calc.get_components()

111 # If any components were merged, create a new interval

112 if k != kprev:

113 # Only keep intervals of non-zero width

114 intervals.append(build_kinterval(kprev, k, hprev,

115 components_prev, cdim_prev))

116 kprev = k

117 hprev = h

118 components_prev = components

119 cdim_prev = cdim

121 # Stop once all components are merged

122 if h == end_state:

123 intervals.append(build_kinterval(k, float("inf"), h,

124 components, cdim))

125 return intervals

128def analyze_kintervals(atoms, method='RDA', merge=True):

129 """Performs a k-interval analysis.

131 In each k-interval the components (connected clusters) are identified.

132 The intervals are sorted according to the scoring parameter, from high

133 to low.

135 Parameters:

137 atoms: ASE atoms object

138 The system to analyze. The periodic boundary conditions determine

139 the maximum achievable component dimensionality, i.e. pbc=[1,1,0] sets

140 a maximum dimensionality of 2.

141 method: string

142 Analysis method to use, either 'RDA' (default option) or 'TSA'.

143 These correspond to the Rank Determination Algorithm of Mounet et al.

144 and the Topological Scaling Algorithm (TSA) of Ashton et al.

145 merge: boolean

146 Decides if k-intervals of the same type (e.g. 01D or 3D) should be

147 merged. Default: true

149 Returns:

151 intervals: list

152 List of KIntervals for each interval identified. A KInterval is a

153 namedtuple with the following field names:

155 score: float

156 Dimensionality score in the range [0, 1]

157 a: float

158 The start of the k-interval

159 b: float

160 The end of the k-interval

161 dimtype: str

162 The dimensionality type

163 h: tuple

164 The histogram of the number of components of each dimensionality.

165 For example, (8, 0, 3, 0) means eight 0D and three 2D components.

166 components: array

167 The component ID of each atom.

168 cdim: dict

169 The component dimensionalities

170 """

171 intervals = build_kintervals(atoms, method)

172 if merge:

173 intervals = merge_intervals(intervals)

175 # Sort intervals by score. Interval width resolves ambiguity when score=0.

176 return sorted(intervals, reverse=True, key=lambda x: (x.score, x.b - x.a))