Révision 10 ase/embed.py
| embed.py (revision 10) | ||
|---|---|---|
| 1 |
import math |
|
| 1 | 2 |
from ase import Atom, Atoms |
| 2 | 3 |
from ase.data import covalent_radii, atomic_numbers |
| 3 | 4 |
|
| ... | ... | |
| 5 | 6 |
|
| 6 | 7 |
class Embed(Atoms): |
| 7 | 8 |
#--- constructor of the Embed class --- |
| 8 |
def __init__(self, system=None, cluster=None): |
|
| 9 |
def __init__(self, system=None, cluster=None, cell_cluster = "Auto"):
|
|
| 9 | 10 |
super(Embed, self).__init__() |
| 10 | 11 |
# define the atom map |
| 11 | 12 |
self.atom_map_sys_cl = [] |
| ... | ... | |
| 19 | 20 |
# define the systems for calculations |
| 20 | 21 |
self.set_system(system) |
| 21 | 22 |
self.set_cluster(cluster) |
| 22 |
# |
|
| 23 |
self.set_cell([10, 10, 10]) |
|
| 23 |
# set the cell of the system |
|
| 24 |
self.set_cell(system.get_cell()) |
|
| 25 |
self.cell_cluster = cell_cluster |
|
| 24 | 26 |
return |
| 25 | 27 |
|
| 26 | 28 |
#--- set the cluster --- |
| ... | ... | |
| 77 | 79 |
self.find_cluster() |
| 78 | 80 |
self.set_linkatoms() |
| 79 | 81 |
print "link atoms found: ", len(self.linkatoms) |
| 82 |
if self.cell_cluster == "System": |
|
| 83 |
self.atoms_cluster.set_cell(self.atoms_system.get_cell()) |
|
| 84 |
elif self.cell_cluster == "Auto": |
|
| 85 |
positions = self.atoms_cluster.get_positions() |
|
| 86 |
#find the biggest dimensions of the cluster in x,y,z direction |
|
| 87 |
l = 0 |
|
| 88 |
for idir in xrange(3): |
|
| 89 |
l = max(l, positions[:, idir].max() - positions[:, idir].min()) |
|
| 90 |
# calculate the box parameters (cluster + min 5 Ang) |
|
| 91 |
l = (math.floor(l/2.5)+1)*2.5 + 5.0 |
|
| 92 |
# build cell |
|
| 93 |
cell = np.zeros((3, 3), float) |
|
| 94 |
# apply cell parameters |
|
| 95 |
for idir in xrange(3): |
|
| 96 |
cell[idir, idir] = l |
|
| 97 |
# set parameters to cluster |
|
| 98 |
self.atoms_cluster.set_cell(cell) |
|
| 99 |
# print information on the screen |
|
| 100 |
print "length of box surrounding the cluster: ", |
|
| 101 |
print l*10, |
|
| 102 |
print "pm" |
|
| 103 |
else: |
|
| 104 |
self.atoms_cluster.set_cell(self.cell_cluster) |
|
| 80 | 105 |
|
| 81 | 106 |
def find_cluster(self): |
| 82 | 107 |
# set tolerance |
| ... | ... | |
| 122 | 147 |
xyzs_sys=[] |
| 123 | 148 |
for atom_sys in self.atoms_system: |
| 124 | 149 |
xyzs_sys.append(atom_sys.get_position()) |
| 125 |
# FIXME: mapping does not work when cluster atoms are outside of the box |
|
| 126 | 150 |
# set the standard link atom |
| 127 | 151 |
if linkAtom is None: |
| 128 | 152 |
linkAtom ='H' |
| ... | ... | |
| 135 | 159 |
bonds = [] |
| 136 | 160 |
# set the 27 cell_vec, starting with the (0,0,0) vector for the unit cell |
| 137 | 161 |
cells_L = [(0.0, 0.0, 0.0)] |
| 138 |
# get the cell vectors |
|
| 139 |
cell = self.get_cell() |
|
| 140 |
# TODO: check 0D, 1D, 2D, 3D |
|
| 141 |
if False: |
|
| 142 |
for ix in xrange(-1, 2): |
|
| 143 |
for iy in xrange(-1, 2): |
|
| 144 |
for iz in xrange(-1, 2): |
|
| 145 |
if ix == 0 and iy == 0 and iz == 0: |
|
| 146 |
continue |
|
| 147 |
|
|
| 148 |
cells_L.append(np.dot([ix, iy, iz], cell)) |
|
| 162 |
# get the cell vectors of the host system and build up a 3 by 3 supercell to search for neighbors in the surrounding |
|
| 163 |
cell = self.atoms_system.get_cell() |
|
| 164 |
if self.atoms_system.get_pbc().any(): |
|
| 165 |
for ix in xrange(-1, 2): |
|
| 166 |
for iy in xrange(-1, 2): |
|
| 167 |
for iz in xrange(-1, 2): |
|
| 168 |
if ix == 0 and iy == 0 and iz == 0: |
|
| 169 |
continue |
|
| 170 |
cells_L.append(np.dot([ix, iy, iz], cell)) |
|
| 149 | 171 |
# save the radius of system atoms |
| 150 | 172 |
rs_sys = [] |
| 151 | 173 |
for atom in self.atoms_system: |
Formats disponibles : Unified diff