r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1import numpy as np

4from ase.geometry import get_distances

7def random_unit_vector(rng):

8 """Random unit vector equally distributed on the sphere

10 Parameter

11 ---------

12 rng: random number generator object

13 """

14 ct = -1 + 2 * rng.random()

15 phi = 2 * np.pi * rng.random()

16 st = np.sqrt(1 - ct**2)

17 return np.array([st * np.cos(phi), st * np.sin(phi), ct])

20def nearest(atoms1, atoms2, cell=None, pbc=None):

21 """Return indices of nearest atoms"""

22 p1 = atoms1.get_positions()

23 p2 = atoms2.get_positions()

24 vd_aac, d2_aa = get_distances(p1, p2, cell, pbc)

25 i1, i2 = np.argwhere(d2_aa == d2_aa.min())[0]

26 return i1, i2, vd_aac[i1, i2]

29def attach(atoms1, atoms2, distance, direction=(1, 0, 0),

30 maxiter=50, accuracy=1e-5):

31 """Attach two structures

33 Parameters

34 ----------

35 atoms1: Atoms

36 cell and pbc of this object are used

37 atoms2: Atoms

38 distance: float

39 minimal distance (Angstrom)

40 direction: unit vector (3 floats)

41 relative direction between center of masses

42 maxiter: int

43 maximal number of iterations to get required distance, default 100

44 accuracy: float

45 required accuracy for minimal distance (Angstrom), default 1e-5

47 Returns

48 -------

49 Joined structure as an atoms object.

50 """

51 atoms = atoms1.copy()

52 atoms2 = atoms2.copy()

54 direction = np.array(direction, dtype=float)

55 direction /= np.linalg.norm(direction)

56 assert len(direction) == 3

57 dist2 = distance**2

59 i1, i2, dv_c = nearest(atoms, atoms2, atoms.cell, atoms.pbc)

61 for i in range(maxiter):

62 dv2 = (dv_c**2).sum()

64 vcost = np.dot(dv_c, direction)

65 a = np.sqrt(max(0, dist2 - dv2 + vcost**2))

66 move = a - vcost

67 if abs(move) < accuracy:

68 atoms += atoms2

69 return atoms

71 # we need to move

72 atoms2.translate(direction * move)

73 i1, i2, dv_c = nearest(atoms, atoms2, atoms.cell, atoms.pbc)

75 raise RuntimeError('attach did not converge')

78def attach_randomly(atoms1, atoms2, distance,

79 rng=np.random):

80 """Randomly attach two structures with a given minimal distance

82 Parameters

83 ----------

84 atoms1: Atoms object

85 atoms2: Atoms object

86 distance: float

87 Required distance

88 rng: random number generator object

89 defaults to np.random.RandomState()

91 Returns

92 -------

93 Joined structure as an atoms object.

94 """

95 atoms2 = atoms2.copy()

96 atoms2.rotate('x', random_unit_vector(rng),

97 center=atoms2.get_center_of_mass())

98 return attach(atoms1, atoms2, distance,

99 direction=random_unit_vector(rng))

103 rng=np.random,

104 comm=world):

105 """Randomly attach two structures with a given minimal distance

106 and ensure that these are distributed.

108 Parameters

109 ----------

110 atoms1: Atoms object

111 atoms2: Atoms object

112 distance: float

113 Required distance

114 rng: random number generator object

115 defaults to np.random.RandomState()

116 comm: communicator to distribute

117 Communicator to distribute the structure, default: world

119 Returns

120 -------

121 Joined structure as an atoms object.

122 """

123 if comm.rank == 0:

124 joined = attach_randomly(atoms1, atoms2, distance, rng)