314 |
314 |
|
315 |
315 |
|
316 |
316 |
def check_collision(slab_molec, slab_num_atoms, min_height, vect, nn_slab=0,
|
317 |
|
nn_molec=0, coll_coeff=1.2):
|
|
317 |
nn_molec=0, coll_coeff=1.2, exclude_atom=False):
|
318 |
318 |
"""Checks whether a slab and a molecule collide or not.
|
319 |
319 |
|
320 |
320 |
@param slab_molec: The system of adsorbate-slab for which to detect if there
|
... | ... | |
328 |
328 |
@param min_height: The minimum height atoms can have to not be considered as
|
329 |
329 |
colliding.
|
330 |
330 |
@param vect: The vector perpendicular to the slab.
|
|
331 |
@param exclude_atom: Whether to exclude the adsorption center in the
|
|
332 |
molecule in the collision detection.
|
331 |
333 |
@return: bool, whether the surface and the molecule collide.
|
332 |
334 |
"""
|
|
335 |
from copy import deepcopy
|
333 |
336 |
from ase.neighborlist import natural_cutoffs, neighbor_list
|
334 |
337 |
|
335 |
338 |
# Check structure overlap by height
|
... | ... | |
342 |
345 |
logger.error(err_msg)
|
343 |
346 |
raise ValueError(err_msg)
|
344 |
347 |
for atom in slab_molec[slab_num_atoms:]:
|
|
348 |
if exclude_atom is not False \
|
|
349 |
and atom.index == exclude_atom:
|
|
350 |
continue
|
345 |
351 |
for i, coord in enumerate(vect):
|
346 |
352 |
if coord == 0:
|
347 |
353 |
continue
|
... | ... | |
350 |
356 |
|
351 |
357 |
# Check structure overlap by sphere collision
|
352 |
358 |
if coll_coeff is not False:
|
353 |
|
slab_molec_cutoffs = natural_cutoffs(slab_molec, mult=coll_coeff)
|
354 |
|
slab_molec_nghbs = len(
|
355 |
|
neighbor_list("i", slab_molec, slab_molec_cutoffs))
|
|
359 |
if exclude_atom is not False:
|
|
360 |
slab_molec_wo_ctr = deepcopy(slab_molec)
|
|
361 |
del slab_molec_wo_ctr[exclude_atom + slab_num_atoms]
|
|
362 |
slab_molec_cutoffs = natural_cutoffs(slab_molec_wo_ctr,
|
|
363 |
mult=coll_coeff)
|
|
364 |
slab_molec_nghbs = len(neighbor_list("i", slab_molec_wo_ctr,
|
|
365 |
slab_molec_cutoffs))
|
|
366 |
else:
|
|
367 |
slab_molec_cutoffs = natural_cutoffs(slab_molec, mult=coll_coeff)
|
|
368 |
slab_molec_nghbs = len(neighbor_list("i", slab_molec,
|
|
369 |
slab_molec_cutoffs))
|
356 |
370 |
if slab_molec_nghbs > nn_slab + nn_molec:
|
357 |
371 |
return True
|
358 |
372 |
|
... | ... | |
361 |
375 |
|
362 |
376 |
def correct_coll(molec, slab, ctr_coord, site_coord, num_pts,
|
363 |
377 |
min_coll_height, norm_vect, slab_nghbs, molec_nghbs,
|
364 |
|
coll_coeff, height=2.5):
|
365 |
|
# TODO Rethink this function
|
|
378 |
coll_coeff, height=2.5, excl_atom=False):
|
|
379 |
# TODO Rename this function
|
366 |
380 |
"""Tries to adsorb a molecule on a slab trying to avoid collisions by doing
|
367 |
381 |
small rotations.
|
368 |
382 |
|
... | ... | |
382 |
396 |
atoms resulting in a distance that two atoms being closer to that is
|
383 |
397 |
considered as atomic collision.
|
384 |
398 |
@param height: Height on which to try adsorption.
|
|
399 |
@param excl_atom: Whether to exclude the adsorption center in the
|
|
400 |
molecule in the collision detection.
|
385 |
401 |
@return collision: bool, whether the structure generated has collisions
|
386 |
402 |
between slab and adsorbate.
|
387 |
403 |
"""
|
... | ... | |
401 |
417 |
norm_vect=norm_vect)
|
402 |
418 |
collision = check_collision(slab_molec, slab_num_atoms, min_coll_height,
|
403 |
419 |
norm_vect, slab_nghbs, molec_nghbs,
|
404 |
|
coll_coeff)
|
|
420 |
coll_coeff, excl_atom)
|
405 |
421 |
num_corr += 1
|
406 |
422 |
return slab_molec, collision
|
407 |
423 |
|
... | ... | |
506 |
522 |
|
507 |
523 |
def ads_euler(orig_molec, slab, ctr_coord, site_coord, num_pts,
|
508 |
524 |
min_coll_height, coll_coeff, norm_vect, slab_nghbs, molec_nghbs,
|
509 |
|
h_donor, h_acceptor, height):
|
|
525 |
h_donor, h_acceptor, height, excl_atom):
|
510 |
526 |
"""Generates adsorbate-surface structures by sampling over Euler angles.
|
511 |
527 |
|
512 |
528 |
This function generates a number of adsorbate-surface structures at
|
... | ... | |
530 |
546 |
@param h_donor: List of atom types or atom numbers that are H-donors.
|
531 |
547 |
@param h_acceptor: List of atom types or atom numbers that are H-acceptors.
|
532 |
548 |
@param height: Height on which to try adsorption.
|
|
549 |
@param excl_atom: Whether to exclude the adsorption center in the
|
|
550 |
molecule in the collision detection.
|
533 |
551 |
@return: list of ase.Atoms object conatining all the orientations of a given
|
534 |
552 |
conformer.
|
535 |
553 |
"""
|
... | ... | |
546 |
564 |
continue
|
547 |
565 |
molec = deepcopy(prealigned_molec)
|
548 |
566 |
molec.euler_rotate(alpha, beta, gamma, center=ctr_coord)
|
549 |
|
slab_molec, collision = correct_coll(molec, slab,
|
550 |
|
ctr_coord, site_coord,
|
551 |
|
num_pts, min_coll_height,
|
552 |
|
norm_vect,
|
|
567 |
slab_molec, collision = correct_coll(molec, slab, ctr_coord,
|
|
568 |
site_coord, num_pts,
|
|
569 |
min_coll_height, norm_vect,
|
553 |
570 |
slab_nghbs, molec_nghbs,
|
554 |
|
coll_coeff, height)
|
|
571 |
coll_coeff, height,
|
|
572 |
excl_atom)
|
555 |
573 |
if not collision and not any([np.allclose(slab_molec.positions,
|
556 |
574 |
conf.positions)
|
557 |
575 |
for conf in slab_ads_list]):
|
... | ... | |
832 |
850 |
|
833 |
851 |
def ads_internal(orig_molec, slab, ctr1_mol, ctr2_mol, ctr3_mol, ctr1_surf,
|
834 |
852 |
ctr2_surf, num_pts, min_coll_height, coll_coeff, norm_vect,
|
835 |
|
slab_nghbs, molec_nghbs, h_donor, h_acceptor, max_hel, height):
|
|
853 |
slab_nghbs, molec_nghbs, h_donor, h_acceptor, max_hel, height,
|
|
854 |
excl_atom):
|
836 |
855 |
"""Generates adsorbate-surface structures by sampling over internal angles.
|
837 |
856 |
|
838 |
857 |
@param orig_molec: ase.Atoms object of the molecule to adsorb (adsorbate).
|
... | ... | |
861 |
880 |
@param max_hel: Maximum value for sampling the helicopter
|
862 |
881 |
(surf-adsorbate dihedral) angle.
|
863 |
882 |
@param height: Height on which to try adsorption.
|
|
883 |
@param excl_atom: Whether to exclude the adsorption center in the
|
|
884 |
molecule in the collision detection.
|
864 |
885 |
@return: list of ase.Atoms object conatining all the orientations of a given
|
865 |
886 |
conformer.
|
866 |
887 |
"""
|
... | ... | |
887 |
908 |
num_pts, min_coll_height,
|
888 |
909 |
norm_vect, slab_nghbs,
|
889 |
910 |
molec_nghbs, coll_coeff,
|
890 |
|
height)
|
|
911 |
height, excl_atom)
|
891 |
912 |
if not collision and \
|
892 |
913 |
not any([np.allclose(slab_molec.positions,
|
893 |
914 |
conf.positions)
|
... | ... | |
913 |
934 |
@param inp_vars: Calculation parameters from input file.
|
914 |
935 |
@return: list of ase.Atoms for the adsorbate-surface structures
|
915 |
936 |
"""
|
|
937 |
from copy import deepcopy
|
916 |
938 |
from ase.neighborlist import natural_cutoffs, neighbor_list
|
917 |
939 |
molec_ctrs = inp_vars['molec_ctrs']
|
918 |
940 |
sites = inp_vars['sites']
|
... | ... | |
921 |
943 |
inp_norm_vect = inp_vars['surf_norm_vect']
|
922 |
944 |
min_coll_height = inp_vars['min_coll_height']
|
923 |
945 |
coll_coeff = inp_vars['collision_threshold']
|
|
946 |
exclude_ads_ctr = inp_vars['exclude_ads_ctr']
|
924 |
947 |
h_donor = inp_vars['h_donor']
|
925 |
948 |
h_acceptor = inp_vars['h_acceptor']
|
926 |
949 |
height = inp_vars['adsorption_height']
|
... | ... | |
941 |
964 |
if inp_vars['pbc_cell'] is not False:
|
942 |
965 |
conf.set_pbc(True)
|
943 |
966 |
conf.set_cell(inp_vars['pbc_cell'])
|
944 |
|
if coll_coeff is not False:
|
945 |
|
conf_cutoffs = natural_cutoffs(conf, mult=coll_coeff)
|
946 |
|
molec_nghbs = len(neighbor_list("i", conf, conf_cutoffs))
|
947 |
|
else:
|
948 |
|
molec_nghbs = 0
|
949 |
967 |
for s, site in enumerate(sites_coords):
|
950 |
968 |
if isinstance(inp_norm_vect, str) and inp_norm_vect == 'auto':
|
951 |
969 |
norm_vect = compute_norm_vect(surf, sites[s],
|
... | ... | |
953 |
971 |
else:
|
954 |
972 |
norm_vect = inp_norm_vect
|
955 |
973 |
for c, ctr in enumerate(molec_ctr_coords):
|
|
974 |
if exclude_ads_ctr and isinstance(molec_ctrs[c], int):
|
|
975 |
exclude_atom = molec_ctrs[c]
|
|
976 |
else:
|
|
977 |
exclude_atom = False
|
|
978 |
if exclude_ads_ctr and not isinstance(molec_ctrs[c], int):
|
|
979 |
logger.warning("'exclude_ads_ctr' only works for atomic"
|
|
980 |
"centers and not for many-atoms "
|
|
981 |
f"barycenters. {molec_ctrs[c]} are not "
|
|
982 |
f"going to be excluded from collison.")
|
|
983 |
if coll_coeff and exclude_atom is not False:
|
|
984 |
conf_wo_ctr = deepcopy(conf)
|
|
985 |
del conf_wo_ctr[exclude_atom]
|
|
986 |
conf_cutoffs = natural_cutoffs(conf_wo_ctr, mult=coll_coeff)
|
|
987 |
molec_nghbs = len(neighbor_list("i", conf_wo_ctr,
|
|
988 |
conf_cutoffs))
|
|
989 |
elif coll_coeff and exclude_atom is False:
|
|
990 |
conf_cutoffs = natural_cutoffs(conf, mult=coll_coeff)
|
|
991 |
molec_nghbs = len(neighbor_list("i", conf, conf_cutoffs))
|
|
992 |
else:
|
|
993 |
molec_nghbs = 0
|
956 |
994 |
if angles == 'euler':
|
957 |
995 |
surf_ads_list.extend(ads_euler(conf, surf, ctr, site,
|
958 |
996 |
num_pts, min_coll_height,
|
959 |
997 |
coll_coeff, norm_vect,
|
960 |
998 |
surf_nghbs, molec_nghbs,
|
961 |
|
h_donor, h_acceptor, height))
|
|
999 |
h_donor, h_acceptor, height,
|
|
1000 |
exclude_atom))
|
962 |
1001 |
elif angles == 'internal':
|
963 |
1002 |
mol_ctr1 = molec_ctrs[c]
|
964 |
1003 |
mol_ctr2 = inp_vars["molec_ctrs2"][c]
|
... | ... | |
973 |
1012 |
coll_coeff, norm_vect,
|
974 |
1013 |
surf_nghbs, molec_nghbs,
|
975 |
1014 |
h_donor, h_acceptor,
|
976 |
|
max_h, height))
|
|
1015 |
max_h, height,
|
|
1016 |
exclude_atom))
|
977 |
1017 |
return surf_ads_list
|
978 |
1018 |
|
979 |
1019 |
|