Chimie Théorique » scripts_chimie4psmn » DockOnSurf

def get_disso_atoms():

    disso_atoms_str = dos_inp.get('Screening', 'disso_atoms', fallback="False")

    disso_atoms = []

    if disso_atoms_str.lower() in turn_false_answers:

        pass

    else:

        for el in disso_atoms_str.split():

            try:

                disso_atoms.append(int(el))

            except ValueError:

                disso_atoms.append(el)

    return disso_atoms

        return_vars['disso_atoms'] = get_disso_atoms()

        slab_molec_nghbs = len(

            neighbor_list("i", slab_molec, slab_molec_cutoffs))

def dissociate_h(slab_molec_orig, h_idx, num_atoms_slab, neigh_cutoff=1):

    # TODO rethink

    """Moves a H atom from a molecule to the surface if the distance is short

    enough.

    @param slab_molec_orig: The ase.Atoms object of the system adsorbate-slab.

    @param h_idx: The index of the hydrogen atom to carry out adsorption of.

    @param num_atoms_slab: The number of atoms of the slab without adsorbate.

    @param neigh_cutoff: half the maximum distance between the surface and the

        H for it to carry out dissociation.

    @return: An ase.Atoms object of the system adsorbate-surface with H

    """

    from copy import deepcopy

    from ase.neighborlist import NeighborList

    slab_molec = deepcopy(slab_molec_orig)

    cutoffs = len(slab_molec) * [neigh_cutoff]

    nl = NeighborList(cutoffs, self_interaction=False, bothways=True)

    nl.update(slab_molec)

    surf_h_vect = np.array([np.infty] * 3)

    for neigh_idx in nl.get_neighbors(h_idx)[0]:

        if neigh_idx < num_atoms_slab:

            dist = np.linalg.norm(slab_molec[neigh_idx].position -

                                  slab_molec[h_idx].position)

            if dist < np.linalg.norm(surf_h_vect):

                surf_h_vect = slab_molec[neigh_idx].position \

                              - slab_molec[h_idx].position

    if np.linalg.norm(surf_h_vect) != np.infty:

        trans_vect = surf_h_vect - surf_h_vect / np.linalg.norm(surf_h_vect)

        slab_molec[h_idx].position = slab_molec[h_idx].position + trans_vect

        return slab_molec

def dissociation(slab_molec, disso_atoms, num_atoms_slab):

    # TODO rethink

    # TODO multiple dissociation

    """Decides which H atoms to dissociate according to a list of atoms.

    Given a list of chemical symbols or atom indices it checks for every atom

    or any of its neighbor if it's a H and calls dissociate_h to try to carry

    out dissociation of that H. For atom indices, it checks both whether

    the atom index or its neighbors are H, for chemical symbols, it only checks

    if there is a neighbor H.

    @param slab_molec: The ase.Atoms object of the system adsorbate-slab.

    @param disso_atoms: The indices or chemical symbols of the atoms

    @param num_atoms_slab:

    @return:

    """

    from ase.neighborlist import natural_cutoffs, NeighborList

    molec = slab_molec[num_atoms_slab:]

    cutoffs = natural_cutoffs(molec)

    nl = NeighborList(cutoffs, self_interaction=False, bothways=True)

    nl.update(molec)

    disso_structs = []

    for el in disso_atoms:

        if isinstance(el, int):

            if molec[el].symbol == 'H':

                disso_struct = dissociate_h(slab_molec, el + num_atoms_slab,

                                            num_atoms_slab)

                if disso_struct is not None:

                    disso_structs.append(disso_struct)

            else:

                for neigh_idx in nl.get_neighbors(el)[0]:

                    if molec[neigh_idx].symbol == 'H':

                        disso_struct = dissociate_h(slab_molec, neigh_idx +

                                                    num_atoms_slab,

                                                    num_atoms_slab)

                        if disso_struct is not None:

                            disso_structs.append(disso_struct)

        else:

            for atom in molec:

                if atom.symbol.lower() == el.lower():

                    for neigh_idx in nl.get_neighbors(atom.index)[0]:

                        if molec[neigh_idx].symbol == 'H':

                            disso_struct = dissociate_h(slab_molec, neigh_idx \

                                                        + num_atoms_slab,

                                                        num_atoms_slab)

                            if disso_struct is not None:

                                disso_structs.append(disso_struct)

    return disso_structs

              min_coll_height, coll_coeff, norm_vect, slab_nghbs, molec_nghbs,

              disso_atoms):

    @param disso_atoms: List of atom types or atom numbers to try to dissociate.

                    slab_ads_list.extend(dissociation(slab_molec, disso_atoms,

                                                      len(slab)))

                 norm_vect, min_coll_height, coll_coeff, disso_atoms):

    @param disso_atoms: List of atom types or atom numbers to try to dissociate.

                                                   surf_nghbs, molec_nghbs,

                                                   disso_atoms))

                                 inp_vars['collision_threshold'],

                                 inp_vars['disso_atoms'])

        self.assertEqual(get_disso_atoms(), True)

        self.assertRaises(ValueError, get_disso_atoms)