226 |
226 |
if atom.position[i] * coord < min_height:
|
227 |
227 |
return True
|
228 |
228 |
return False
|
229 |
|
else: # TODO Implement PBC for sphere collision
|
|
229 |
else:
|
230 |
230 |
slab_molec_cutoffs = natural_cutoffs(slab_molec, mult=coll_coeff)
|
231 |
231 |
slab_molec_nghbs = len(
|
232 |
232 |
neighbor_list("i", slab_molec, slab_molec_cutoffs))
|
... | ... | |
310 |
310 |
if atom.symbol.lower() == el.lower():
|
311 |
311 |
for neigh_idx in nl.get_neighbors(atom.index)[0]:
|
312 |
312 |
if molec[neigh_idx].symbol == 'H':
|
313 |
|
disso_struct = dissociate_h(slab_molec, neigh_idx \
|
|
313 |
disso_struct = dissociate_h(slab_molec, neigh_idx
|
314 |
314 |
+ num_atoms_slab,
|
315 |
315 |
num_atoms_slab)
|
316 |
316 |
if disso_struct is not None:
|
... | ... | |
346 |
346 |
"""
|
347 |
347 |
from copy import deepcopy
|
348 |
348 |
slab_num_atoms = len(slab)
|
|
349 |
slab_molec = []
|
349 |
350 |
collision = True
|
350 |
351 |
max_corr = 6 # Should be an even number
|
351 |
352 |
d_angle = 180 / ((max_corr / 2.0) * num_pts)
|
... | ... | |
372 |
373 |
This function generates a number of adsorbate-surface structures at
|
373 |
374 |
different orientations of the adsorbate sampled at multiple Euler (zxz)
|
374 |
375 |
angles.
|
375 |
|
@param orig_molec: ase.Atoms object of the molecule to adsorb
|
376 |
|
@param slab: ase.Atoms object of the surface on which to adsorb the molecule
|
|
376 |
@param orig_molec: ase.Atoms object of the molecule to adsorb.
|
|
377 |
@param slab: ase.Atoms object of the surface on which to adsorb the
|
|
378 |
molecule.
|
377 |
379 |
@param ctr_coord: The coordinates of the molecule to use as adsorption
|
378 |
380 |
center.
|
379 |
381 |
@param site_coord: The coordinates of the surface on which to adsorb the
|
380 |
|
molecule
|
|
382 |
molecule.
|
381 |
383 |
@param num_pts: Number on which to sample Euler angles.
|
382 |
|
@param min_coll_height: The lowermost height for which to detect a collision
|
|
384 |
@param min_coll_height: The lowest height for which to detect a collision.
|
383 |
385 |
@param coll_coeff: The coefficient that multiplies the covalent radius of
|
384 |
386 |
atoms resulting in a distance that two atoms being closer to that is
|
385 |
387 |
considered as atomic collision.
|
... | ... | |
388 |
390 |
@param molec_nghbs: Number of neighbors in the molecule.
|
389 |
391 |
@param disso_atoms: List of atom types or atom numbers to try to dissociate.
|
390 |
392 |
@return: list of ase.Atoms object conatining all the orientations of a given
|
391 |
|
conformer
|
|
393 |
conformer.
|
392 |
394 |
"""
|
393 |
395 |
from copy import deepcopy
|
394 |
396 |
slab_ads_list = []
|
... | ... | |
406 |
408 |
norm_vect,
|
407 |
409 |
slab_nghbs, molec_nghbs,
|
408 |
410 |
coll_coeff)
|
409 |
|
if not collision and \
|
410 |
|
not any([(slab_molec.positions == conf.positions).all()
|
411 |
|
for conf in slab_ads_list]):
|
|
411 |
if not collision and not any([np.allclose(slab_molec.positions,
|
|
412 |
conf.positions)
|
|
413 |
for conf in slab_ads_list]):
|
412 |
414 |
slab_ads_list.append(slab_molec)
|
413 |
415 |
slab_ads_list.extend(dissociation(slab_molec, disso_atoms,
|
414 |
416 |
len(slab)))
|
... | ... | |
619 |
621 |
# Check bond rotation is unmodified
|
620 |
622 |
bond_angle = get_vect_angle(-vect_inter, vect_mol)
|
621 |
623 |
if np.isclose((dihedral_angle - dihedral_angle_target + 90) % 180 - 90,
|
622 |
|
0, atol=1e-3) and np.isclose(
|
623 |
|
(bond_angle - bond_angle_target + 90) %
|
624 |
|
180 - 90, 0, atol=1e-5) and np.allclose(get_atom_coords(molecule,
|
|
624 |
0, atol=1e-3) \
|
|
625 |
and np.isclose((bond_angle - bond_angle_target + 90) %
|
|
626 |
180 - 90, 0, atol=1e-5) \
|
|
627 |
and np.allclose(get_atom_coords(molecule,
|
625 |
628 |
ctr1_mol)
|
626 |
629 |
- get_atom_coords(surf,
|
627 |
630 |
ctr1_surf),
|
... | ... | |
645 |
648 |
# vect_mol)
|
646 |
649 |
vect_mol_inner = np.dot(vect_mol, vect_mol)
|
647 |
650 |
bond_inter_reject = -vect_inter - np.dot(-vect_inter, vect_mol) / \
|
648 |
|
vect_mol_inner * vect_mol
|
|
651 |
vect_mol_inner * vect_mol
|
649 |
652 |
bond2_mol_reject = vect2_mol - np.dot(vect2_mol, vect_mol) / \
|
650 |
|
vect_mol_inner * vect_mol
|
|
653 |
vect_mol_inner * vect_mol
|
651 |
654 |
dihedral_angle_ini = get_vect_angle(bond_inter_reject,
|
652 |
655 |
bond2_mol_reject, vect_mol)
|
653 |
656 |
|
... | ... | |
657 |
660 |
|
658 |
661 |
# Update molecular bonds
|
659 |
662 |
vect_mol = get_atom_coords(molecule, ctr2_mol) \
|
660 |
|
- get_atom_coords(molecule, ctr1_mol)
|
|
663 |
- get_atom_coords(molecule, ctr1_mol)
|
661 |
664 |
vect2_mol = get_atom_coords(molecule, ctr3_mol) \
|
662 |
|
- get_atom_coords(molecule, ctr2_mol)
|
|
665 |
- get_atom_coords(molecule, ctr2_mol)
|
663 |
666 |
|
664 |
667 |
# Check if rotation was successful
|
665 |
668 |
# Check adsorbate dihedral rotation
|
666 |
669 |
vect_mol_inner = np.dot(vect_mol, vect_mol)
|
667 |
670 |
bond2_mol_reject = vect2_mol - np.dot(vect2_mol, vect_mol) / \
|
668 |
|
vect_mol_inner * vect_mol
|
|
671 |
vect_mol_inner * vect_mol
|
669 |
672 |
mol_dihedral_angle = get_vect_angle(bond_inter_reject,
|
670 |
673 |
bond2_mol_reject, vect_mol)
|
671 |
674 |
# Check dihedral rotation
|
672 |
675 |
vect_inter_inner = np.dot(vect_inter, vect_inter)
|
673 |
676 |
vect_surf_reject = vect_surf - np.dot(vect_surf, vect_inter) / \
|
674 |
|
vect_inter_inner * vect_inter
|
|
677 |
vect_inter_inner * vect_inter
|
675 |
678 |
vect_mol_reject = vect_mol - np.dot(vect_mol, vect_inter) / \
|
676 |
|
vect_inter_inner * vect_inter
|
|
679 |
vect_inter_inner * vect_inter
|
677 |
680 |
dihedral_angle = get_vect_angle(vect_surf_reject, vect_mol_reject,
|
678 |
681 |
vect_inter)
|
679 |
682 |
# Check bond rotation is unmodified
|
... | ... | |
702 |
705 |
def ads_chemcat(orig_molec, slab, ctr1_mol, ctr2_mol, ctr3_mol, ctr1_surf,
|
703 |
706 |
ctr2_surf, num_pts, min_coll_height, coll_coeff, norm_vect,
|
704 |
707 |
slab_nghbs, molec_nghbs, disso_atoms, max_h):
|
|
708 |
"""Generates adsorbate-surface structures by sampling over chemcat angles.
|
|
709 |
|
|
710 |
@param orig_molec: ase.Atoms object of the molecule to adsorb (adsorbate).
|
|
711 |
@param slab: ase.Atoms object of the surface on which to adsorb the molecule
|
|
712 |
@param ctr1_mol: The index/es of the center in the adsorbate to use as
|
|
713 |
adsorption center.
|
|
714 |
@param ctr2_mol: The index/es of the center in the adsorbate to use for the
|
|
715 |
definition of the surf-adsorbate angle, surf-adsorbate dihedral angle
|
|
716 |
and adsorbate dihedral angle.
|
|
717 |
@param ctr3_mol: The index/es of the center in the adsorbate to use for the
|
|
718 |
definition of the adsorbate dihedral angle.
|
|
719 |
@param ctr1_surf: The index/es of the center in the surface to use as
|
|
720 |
adsorption center.
|
|
721 |
@param ctr2_surf: The index/es of the center in the surface to use for the
|
|
722 |
definition of the surf-adsorbate dihedral angle.
|
|
723 |
@param num_pts: Number on which to sample Euler angles.
|
|
724 |
@param min_coll_height: The lowest height for which to detect a collision
|
|
725 |
@param coll_coeff: The coefficient that multiplies the covalent radius of
|
|
726 |
atoms resulting in a distance that two atoms being closer to that is
|
|
727 |
considered as atomic collision.
|
|
728 |
@param norm_vect: The vector perpendicular to the surface.
|
|
729 |
@param slab_nghbs: Number of neigbors in the surface.
|
|
730 |
@param molec_nghbs: Number of neighbors in the molecule.
|
|
731 |
@param disso_atoms: List of atom types or atom numbers to try to dissociate.
|
|
732 |
@param max_h: Maximum value for sampling the helicopter
|
|
733 |
(surf-adsorbate dihedral) angle.
|
|
734 |
@return: list of ase.Atoms object conatining all the orientations of a given
|
|
735 |
conformer.
|
|
736 |
"""
|
705 |
737 |
from copy import deepcopy
|
706 |
738 |
slab_ads_list = []
|
707 |
739 |
# Rotation over bond angle
|
708 |
740 |
for alpha in np.arange(90, 180, 90 / min(num_pts, 2)):
|
709 |
741 |
# Rotation over surf-adsorbate dihedral angle
|
710 |
|
for beta in np.arange(0, max_h, max_h / num_pts): # Maybe 360
|
|
742 |
for beta in np.arange(0, max_h, max_h / num_pts):
|
711 |
743 |
# Rotation over adsorbate bond dihedral angle
|
712 |
744 |
for gamma in np.arange(0, 360, 360 / num_pts):
|
713 |
745 |
new_molec = deepcopy(orig_molec)
|
... | ... | |
723 |
755 |
molec_nghbs, coll_coeff,
|
724 |
756 |
2.5)
|
725 |
757 |
if not collision and \
|
726 |
|
not any([(slab_molec.positions == conf.positions).all()
|
|
758 |
not any([np.allclose(slab_molec.positions,
|
|
759 |
conf.positions)
|
727 |
760 |
for conf in slab_ads_list]):
|
728 |
761 |
slab_ads_list.append(slab_molec)
|
729 |
762 |
slab_ads_list.extend(dissociation(slab_molec, disso_atoms,
|