Genetic algorithm

Optimize adsorbate coverage pattern

Adsorbate procreation operators that adds an adsorbate to the surface of a particle or given structure.

class acat.ga.adsorbate_operators.AdsorbateOperator(adsorbate_species, num_muts=1)[source]

Bases: ase.ga.offspring_creator.OffspringCreator

Base class for all operators that add, move or remove adsorbates.

Don’t use this operator directly!

classmethod initialize_individual(parent, indi=None)[source]

Initializes a new individual that inherits some parameters from the parent, and initializes the info dictionary. If the new individual already has more structure it can be supplied in the parameter indi.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

add_adsorbate(atoms, hetero_site_list, heights, adsorbate_species=None, min_adsorbate_distance=2.0, tilt_angle=0.0)[source]

Adds the adsorbate in self.adsorbate to the supplied atoms object at the first free site in the specified site_list. A site is free if no other adsorbates can be found in a sphere of radius min_adsorbate_distance around the chosen site.

Parameters
  • atoms (ase.Atoms object) – The atoms object that the adsorbate will be added to.

  • hetero_site_list (list) – A list of dictionaries, each dictionary contains site information given by acat.adsorbate_coverage module.

  • heights (dict) – A dictionary that contains the adsorbate height for each site type.

  • adsorbate_species (str or list of strs, default None) – One or a list of adsorbate species to be added to the surface. Use self.adsorbate_species if not specified.

  • min_adsorbate_distance (float, default 2.) – The radius of the sphere inside which no other adsorbates should be found.

  • tilt_angle (float, default 0.) – Tilt the adsorbate with an angle (in degress) relative to the surface normal.

remove_adsorbate(atoms, hetero_site_list, return_site_index=False)[source]

Removes an adsorbate from the atoms object at the first occupied site in hetero_site_list. If no adsorbates can be found, one will be added instead.

Parameters
  • atoms (ase.Atoms object) – The atoms object that the adsorbate will be added to

  • hetero_site_list (list) – A list of dictionaries, each dictionary contains site information given by acat.adsorbate_coverage module.

  • return_site_index (bool, default False) – Whether to return the site index of the hetero_site_list instead of the site. Useful for moving or replacing adsorbate.

get_adsorbate_indices(atoms, position)[source]

Returns the indices of the adsorbate at the supplied position.

class acat.ga.adsorbate_operators.AddAdsorbate(adsorbate_species, heights={'3fold': 1.3, '4fold': 1.3, '5fold': 1.5, '6fold': 0.0, 'bridge': 1.5, 'fcc': 1.3, 'hcp': 1.3, 'longbridge': 1.5, 'ontop': 1.8, 'shortbridge': 1.5}, min_adsorbate_distance=2.0, surface=None, allow_6fold=False, composition_effect=True, site_preference=None, surface_preference=None, tilt_angle=None, num_muts=1, dmax=2.5)[source]

Bases: acat.ga.adsorbate_operators.AdsorbateOperator

Use this operator to add adsorbates to the surface. The surface is allowed to change during the algorithm run.

There is no limit of adsorbate species. You can either provide one species or a list of species.

Site and surface preference can be supplied. If both are supplied site will be considered first.

Supplying a tilt angle will tilt the adsorbate with an angle relative to the standard perpendicular to the surface.

The operator is generalized for both periodic and non-periodic systems (distinguished by atoms.pbc).

Parameters
  • adsorbate_species (str or list of strs) – One or a list of adsorbate species to be added to the surface.

  • heights (dict, default acat.settings.site_heights) – A dictionary that contains the adsorbate height for each site type. Use the default height settings if the height for a site type is not specified.

  • min_adsorbate_distance (float, default 2.) – The radius of the sphere inside which no other adsorbates should be found.

  • surface (str, default None) – The surface type (crystal structure + Miller indices). Only required if the structure is a periodic surface slab.

  • allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.

  • composition_effect (bool, default True) – Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition=False for monometallics. Since GA is commonly used for bimetallics, the default is set to True.

  • site_preference (str, defualt None) – The site type that has higher priority to attach adsorbates.

  • surface_preference (str, default None) – The surface type that has higher priority to attach adsorbates. Please only use this for nanoparticles.

  • tilt_angle (float, default 0.) – Tilt the adsorbate with an angle (in degress) relative to the surface normal.

  • num_muts (int, default 1) – The number of times to perform this operation.

  • dmax (float, default 2.5) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.

get_new_individual(parents)[source]

Returns the new individual as an atoms object.

class acat.ga.adsorbate_operators.RemoveAdsorbate(adsorbate_species, surface=None, allow_6fold=False, composition_effect=True, site_preference=None, surface_preference=None, num_muts=1, dmax=2.5)[source]

Bases: acat.ga.adsorbate_operators.AdsorbateOperator

This operator removes an adsorbate from the surface. It works exactly (but doing the opposite) as the AddAdsorbate operator.

The operator is generalized for both periodic and non-periodic systems (distinguished by atoms.pbc).

Parameters
  • adsorbate_species (str or list of strs) – One or a list of adsorbate species to be removed from the surface.

  • surface (str, default None) – The surface type (crystal structure + Miller indices). Only required if the structure is a periodic surface slab.

  • allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.

  • composition_effect (bool, default True) – Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition=False for monometallics. Since GA is commonly used for bimetallics, the default is set to True.

  • site_preference (str, defualt None) – The site type that has higher priority to attach adsorbates.

  • surface_preference (str, default None) – The surface type that has higher priority to attach adsorbates. Please only use this for nanoparticles.

  • tilt_angle (float, default 0.) – Tilt the adsorbate with an angle (in degress) relative to the surface normal.

  • num_muts (int, default 1) – The number of times to perform this operation.

  • dmax (float, default 2.5) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

class acat.ga.adsorbate_operators.MoveAdsorbate(adsorbate_species, heights={'3fold': 1.3, '4fold': 1.3, '5fold': 1.5, '6fold': 0.0, 'bridge': 1.5, 'fcc': 1.3, 'hcp': 1.3, 'longbridge': 1.5, 'ontop': 1.8, 'shortbridge': 1.5}, min_adsorbate_distance=2.0, surface=None, allow_6fold=False, composition_effect=True, site_preference_from=None, surface_preference_from=None, site_preference_to=None, surface_preference_to=None, tilt_angle=None, num_muts=1, dmax=2.5)[source]

Bases: acat.ga.adsorbate_operators.AdsorbateOperator

This operator removes an adsorbate from the surface and adds it again to a different site, i.e. effectively moving the adsorbate.

The operator is generalized for both periodic and non-periodic systems (distinguished by atoms.pbc).

Parameters
  • adsorbate_species (str or list of strs) – One or a list of adsorbate species to be added to the surface.

  • heights (dict, default acat.settings.site_heights) – A dictionary that contains the adsorbate height for each site type. Use the default height settings if the height for a site type is not specified.

  • min_adsorbate_distance (float, default 2.) – The radius of the sphere inside which no other adsorbates should be found.

  • surface (str, default None) – The surface type (crystal structure + Miller indices). Only required if the structure is a periodic surface slab.

  • allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.

  • composition_effect (bool, default True) – Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition=False for monometallics. Since GA is commonly used for bimetallics, the default is set to True.

  • site_preference_from (str, defualt None) – The site type that has higher priority to detach adsorbates.

  • surface_preference_from (str, default None) – The surface type that has higher priority to detach adsorbates. Please only use this for nanoparticles.

  • site_preference_to (str, defualt None) – The site type that has higher priority to attach adsorbates.

  • surface_preference_to (str, default None) – The surface type that has higher priority to attach adsorbates. Please only use this for nanoparticles.

  • tilt_angle (float, default 0.) – Tilt the adsorbate with an angle (in degress) relative to the surface normal.

  • num_muts (int, default 1) – The number of times to perform this operation.

  • dmax (float, default 2.5) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

class acat.ga.adsorbate_operators.ReplaceAdsorbate(adsorbate_species, heights={'3fold': 1.3, '4fold': 1.3, '5fold': 1.5, '6fold': 0.0, 'bridge': 1.5, 'fcc': 1.3, 'hcp': 1.3, 'longbridge': 1.5, 'ontop': 1.8, 'shortbridge': 1.5}, min_adsorbate_distance=2.0, surface=None, allow_6fold=False, composition_effect=True, site_preference_from=None, surface_preference_from=None, tilt_angle=None, num_muts=1, dmax=2.5)[source]

Bases: acat.ga.adsorbate_operators.AdsorbateOperator

This operator removes an adsorbate from the surface and adds another species to the same site, i.e. effectively replacing the adsorbate.

The operator is generalized for both periodic and non-periodic systems (distinguished by atoms.pbc).

Parameters
  • adsorbate_species (str or list of strs) – One or a list of adsorbate species to be added to the surface.

  • heights (dict, default acat.settings.site_heights) – A dictionary that contains the adsorbate height for each site type. Use the default height settings if the height for a site type is not specified.

  • min_adsorbate_distance (float, default 2.) – The radius of the sphere inside which no other adsorbates should be found.

  • surface (str, default None) – The surface type (crystal structure + Miller indices). Only required if the structure is a periodic surface slab.

  • allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.

  • composition_effect (bool, default True) – Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition=False for monometallics. Since GA is commonly used for bimetallics, the default is set to True.

  • site_preference_from (str, defualt None) – The site type that has higher priority to replace adsorbates.

  • surface_preference_from (str, default None) – The surface type that has higher priority to replace adsorbates. Please only use this for nanoparticles.

  • tilt_angle (float, default 0.) – Tilt the adsorbate with an angle (in degress) relative to the surface normal.

  • num_muts (int, default 1) – The number of times to perform this operation.

  • dmax (float, default 2.5) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

class acat.ga.adsorbate_operators.CutSpliceCrossoverWithAdsorbates(adsorbate_species, blmin, keep_composition=True, fix_coverage=False, min_adsorbate_distance=2.0, allow_6fold=False, composition_effect=True, rotate_vectors=None, rotate_angles=None, dmax=2.5)[source]

Bases: acat.ga.adsorbate_operators.AdsorbateOperator

Crossover that cuts two particles with adsorbates through a plane in space and merges two halfes from different particles together.

It keeps the correct composition by randomly assigning elements in the new particle. If some of the atoms in the two particle halves are too close, the halves are moved away from each other perpendicular to the cutting plane.

The complexity of crossover with adsorbates makes this operator not robust enough. The adsorption site identification will fail once the nanoparticle shape becomes too irregular after crossover.

Parameters
  • adsorbate_species (str or list of strs) – One or a list of adsorbate species to be added to the surface.

  • blmin (dict) – Dictionary of minimum distance between atomic numbers. e.g. {(28,29): 1.5}

  • keep_composition (bool, default True) – Should the composition be the same as in the parents.

  • fix_coverage (bool, default False) – Should the adsorbate coverage be the same as in the parents.

  • min_adsorbate_distance (float, default 2.) – The radius of the sphere inside which no other adsorbates should be found.

  • allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.

  • composition_effect (bool, default True) – Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition=False for monometallics. Since GA is commonly used for bimetallics, the default is set to True.

  • rotate_vectors (list, default None) – A list of vectors that the part of the structure that is cut is able to rotate around, the size of rotation is set in rotate_angles. Default None meaning no rotation is performed.

  • rotate_angles (list, default None) – A list of angles that the structure cut can be rotated. The vector being rotated around is set in rotate_vectors. Default None meaning no rotation is performed.

  • dmax (float, default 2.5) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

get_vectors_below_min_dist(atoms)[source]

Generator function that returns each vector (between atoms) that is shorter than the minimum distance for those atom types (set during the initialization in blmin).

Comparator objects relevant to particles with adsorbates.

acat.ga.adsorbate_comparators.count_ads(atoms, adsorbate)[source]

Very naive implementation only taking into account the symbols. atoms and adsorbate should both be supplied as Atoms objects.

class acat.ga.adsorbate_comparators.AdsorbateCountComparator(adsorbate)[source]

Bases: object

Compares the number of adsorbates on the particles and returns True if the numbers are the same, False otherwise.

Parameters:

adsorbate: list or string a supplied list of adsorbates or a string if only one adsorbate is possible

looks_like(a1, a2)[source]

Does the actual comparison.

class acat.ga.adsorbate_comparators.AdsorptionSitesComparator(min_diff_adsorption_sites=2)[source]

Bases: object

Compares the metal atoms in the adsorption sites and returns True if less than min_diff_adsorption_sites of the sites with adsorbates consist of different atoms.

Ex: a1.info[‘data’][‘adsorbates_site_atoms’] = [(‘Cu’,’Ni’),(‘Cu’,’Ni’),(‘Ni’),(‘Ni’)]

a2.info[‘data’][‘adsorbates_site_atoms’] = [(‘Cu’,’Ni’),(‘Ni’,’Ni’, ‘Ni’),(‘Ni’),(‘Ni’)]

will have a difference of 2: (2*(‘Cu’,’Ni’)-1*(‘Cu’,’Ni’)=1, 1*(‘Ni’,’Ni’,’Ni’)=1, 2*(‘Ni’)-2*(‘Ni’)=0)

looks_like(a1, a2)[source]
class acat.ga.adsorbate_comparators.AdsorptionMetalsComparator(same_adsorption_number)[source]

Bases: object

Compares the number of adsorbate-metal bonds and returns True if the number for a1 and a2 differs by less than the supplied parameter same_adsorption_number

Ex: a1.info[‘data’][‘adsorbates_bound_to’] = {‘Cu’:1, ‘Ni’:3} a2.info[‘data’][‘adsorbates_bound_to’] = {‘Cu’:.5, ‘Ni’:3.5} will have a difference of .5 in both elements:

looks_like(a1, a2)[source]

Example

All the adsorbate operators can be easily used with other ASE operators. AddAdsorbate, RemoveAdsorbate, MoveAdsorbate and ReplaceAdsorbate operators can be used for both non-periodic nanoparticles and periodic surface slabs. CutSpliceCrossoverWithAdsorbates operator only works for nanoparticles, and it is not recommonded as it is not stable yet.

As an example we will optimize both the adsorbate coverage pattern and the catalyst chemical ordering of a Ni110Pt37 icosahedral nanoalloy with adsorbate species of H, C, O, OH, CO, CH, CH2 and CH3 using the EMT calculator.

The script for a parallel genetic algorithm looks as follows:

from acat.settings import adsorbate_elements
from acat.adsorbate_coverage import ClusterAdsorbateCoverage
from acat.build.ordering import RandomOrderingGenerator as ROG
from acat.build.pattern import random_coverage_pattern
from acat.ga.adsorbate_operators import AddAdsorbate, RemoveAdsorbate
from acat.ga.adsorbate_operators import MoveAdsorbate, ReplaceAdsorbate
from acat.ga.adsorbate_operators import CutSpliceCrossoverWithAdsorbates
from ase.ga.particle_mutations import RandomPermutation, COM2surfPermutation
from ase.ga.particle_mutations import Rich2poorPermutation, Poor2richPermutation
from acat.ga.adsorbate_comparators import AdsorptionSitesComparator
from ase.ga.particle_comparator import NNMatComparator
from ase.ga.standard_comparators import SequentialComparator
from ase.ga.offspring_creator import OperationSelector
from ase.ga.population import Population, RankFitnessPopulation
from ase.ga.convergence import GenerationRepetitionConvergence
from ase.ga.utilities import closest_distances_generator
from ase.ga.data import DataConnection, PrepareDB
from ase.io import read, write
from ase.cluster import Icosahedron
from ase.calculators.emt import EMT
from ase.optimize import BFGS
from collections import defaultdict
from random import choices, uniform
from multiprocessing import Pool
import os

# Define population
# Recommand to choose a number that is a multiple of the number of cpu
pop_size = 50

# Generate 50 icosahedral Ni110Pt37 nanoparticles with random orderings
particle = Icosahedron('Ni', noshells=4)
particle.center(vacuum=5.)
rog = ROG(particle, elements=['Ni', 'Pt'],
          composition={'Ni': 0.75, 'Pt': 0.25},
          trajectory='starting_generation.traj')
rog.run(n_gen=pop_size)

# Generate random coverage on each nanoparticle
species=['H', 'C', 'O', 'OH', 'CO', 'CH', 'CH2', 'CH3']
images = read('starting_generation.traj', index=':')
patterns = []
for atoms in images:
    dmin = uniform(3.5, 8.5)
    pattern = random_coverage_pattern(atoms, adsorbate_species=species,
                                      min_adsorbate_distance=dmin)
    patterns.append(pattern)

# Instantiate the db
db_name = 'ridge_Ni110Pt37_ads.db'

db = PrepareDB(db_name, cell=particle.cell, population_size=pop_size)

for atoms in patterns:
    if 'data' not in atoms.info:
        atoms.info['data'] = {}
    db.add_unrelaxed_candidate(atoms, data=atoms.info['data'])

# Connect to the db
db = DataConnection(db_name)

# Define operators
soclist = ([1, 1, 2, 1, 1, 1, 1],
           [Rich2poorPermutation(elements=['Ni', 'Pt'], num_muts=5),
            Poor2richPermutation(elements=['Ni', 'Pt'], num_muts=5),
            RandomPermutation(num_muts=5),
            AddAdsorbate(species, num_muts=5),
            RemoveAdsorbate(species, num_muts=5),
            MoveAdsorbate(species, num_muts=5),
            ReplaceAdsorbate(species, num_muts=5),])

op_selector = OperationSelector(*soclist)

# Define comparators
comp = SequentialComparator([AdsorptionSitesComparator(10),
                             NNMatComparator(0.2, ['Ni', 'Pt'])],
                            [0.5, 0.5])

def get_ads(atoms):
    """Returns a list of adsorbate names and corresponding indices."""

    if 'data' not in atoms.info:
        atoms.info['data'] = {}
    if 'adsorbates' in atoms.info['data']:
        adsorbates = atoms.info['data']['adsorbates']
    else:
        cac = ClusterAdsorbateCoverage(atoms)
        adsorbates = cac.get_adsorbates()

    return adsorbates

def vf(atoms):
    """Returns the descriptor that distinguishes candidates in the
    niched population."""

    return len(get_ads(atoms))

# Give fittest candidates at different coverages equal fitness
pop = RankFitnessPopulation(data_connection=db,
                            population_size=pop_size,
                            comparator=comp,
                            variable_function=vf,
                            exp_function=True,
                            logfile='log.txt')

# Normal fitness ranking regardless of coverage
#pop = Population(data_connection=db,
#                 population_size=pop_size,
#                 comparator=comp,
#                 logfile='log.txt')

# Set convergence criteria
cc = GenerationRepetitionConvergence(pop, 5)

# Calculate chemical potentials
chem_pots = {'CH4': -24.039, 'H2O': -14.169, 'H2': -6.989}

# Define the relax function
def relax(atoms, single_point=False):
    atoms.center(vacuum=5.)
    atoms.calc = EMT()
    if not single_point:
        opt = BFGS(atoms, logfile=None)
        opt.run(fmax=0.1)

    Epot = atoms.get_potential_energy()
    num_H = len([s for s in atoms.symbols if s == 'H'])
    num_C = len([s for s in atoms.symbols if s == 'C'])
    num_O = len([s for s in atoms.symbols if s == 'O'])
    mutot = num_C * chem_pots['CH4'] + num_O * chem_pots['H2O'] + (
            num_H - 4 * num_C - 2 * num_O) * chem_pots['H2'] / 2
    f = Epot - mutot

    atoms.info['key_value_pairs']['raw_score'] = -f
    atoms.info['key_value_pairs']['potential_energy'] = Epot

    return atoms

# Relax starting generation
def relax_an_unrelaxed_candidate(atoms):
    if 'data' not in atoms.info:
        atoms.info['data'] = {}
    nncomp = atoms.get_chemical_formula(mode='hill')
    print('Relaxing ' + nncomp)
    relax(atoms, single_point=True) # Single point for simplification
    db.add_relaxed_step(atoms)

# Create a multiprocessing Pool
pool = Pool(os.cpu_count())
# Perform relaxations in parallel. Especially
# useful when running GA on large nanoparticles
pool.map(relax_an_unrelaxed_candidate, db.get_all_unrelaxed_candidates())
pop.update()

# Number of generations
num_gens = 1000

# Below is the iterative part of the algorithm
gen_num = db.get_generation_number()
for i in range(num_gens):
    # Check if converged
    if cc.converged():
        print('Converged')
        break
    print('Creating and evaluating generation {0}'.format(gen_num + i))

    def procreation(x):
        # Select an operator and use it
        op = op_selector.get_operator()
        # Select parents for a new candidate
        p1, p2 = pop.get_two_candidates()
        parents = [p1, p2]
        # Pure or bare nanoparticles are not considered
        if len(set(p1.numbers)) < 3:
            continue
        offspring, desc = op.get_new_individual(parents)
        # An operator could return None if an offspring cannot be formed
        # by the chosen parents
        if offspring is None:
            continue
        nncomp = offspring.get_chemical_formula(mode='hill')
        print('Relaxing ' + nncomp)
        if 'data' not in offspring.info:
            offspring.info['data'] = {}
        relax(offspring, single_point=True) # Single point for simplification
        db.add_relaxed_candidate(offspring)

    # Create a multiprocessing Pool
    pool = Pool(os.cpu_count())
    # Perform procreations in parallel. Especially useful when
    # using adsorbate operators which requires site identification
    pool.map(procreation, range(pop_size))

    # update the population to allow new candidates to enter
    pop.update()

Symmetry-constrained genetic algorithm for nanoalloys

Procreation operators meant to be used in symmetry-constrained genetic algorithm (SCGA) for alloy particles.

class acat.ga.symmetric_particle_operators.Mutation(num_muts=1)[source]

Bases: ase.ga.offspring_creator.OffspringCreator

Base class for all particle mutation type operators. Do not call this class directly.

class acat.ga.symmetric_particle_operators.SymmetricSubstitute(shells, elements=None, num_muts=1)[source]

Bases: acat.ga.symmetric_particle_operators.Mutation

Substitute all the atoms in a shell with a different metal element. The elemental composition cannot be fixed.

Parameters
  • shells (list of lists) – The atom indices in each shell divided by symmetry. Can be obtained by acat.build.orderings.SymmetricOrderingGenerator.

  • elements (list of strs, default None) – The metal elements of the nanoalloy. Only take into account the elements specified in this list. Default is to take all elements into account.

  • num_muts (int, default 1) – The number of times to perform this operation.

substitute(atoms)[source]

Does the actual substitution

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

class acat.ga.symmetric_particle_operators.SymmetricPermutation(shells, elements=None, keep_composition=False, num_muts=1)[source]

Bases: acat.ga.symmetric_particle_operators.Mutation

Permutes the elements in two random shells. The elemental composition can be fixed.

Parameters
  • shells (list of lists) – The atom indices in each shell divided by symmetry. Can be obtained by acat.build.orderings.SymmetricOrderingGenerator.

  • elements (list of strs, default None) – The metal elements of the nanoalloy. Only take into account the elements specified in this list. Default is to take all elements into account.

  • keep_composition (bool, defulat False) – Whether the elemental composition should be the same as in the parents.

  • num_muts (int, default 1) – The number of times to perform this operation.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

classmethod mutate(atoms, shells, elements=None, keep_composition=False)[source]

Do the actual permutation.

class acat.ga.symmetric_particle_operators.Crossover[source]

Bases: ase.ga.offspring_creator.OffspringCreator

Base class for all particle crossovers. Do not call this class directly.

class acat.ga.symmetric_particle_operators.SymmetricCrossover(shells, elements=None, keep_composition=False)[source]

Bases: acat.ga.symmetric_particle_operators.Crossover

Merge the elemental distributions in two half shells from different particles together. The elemental composition can be fixed.

Parameters
  • shells (list of lists) – The atom indices in each shell divided by symmetry. Can be obtained by acat.build.orderings.SymmetricOrderingGenerator.

  • elements (list of strs, default None) – The metal elements of the nanoalloy. Only take into account the elements specified in this list. Default is to take all elements into account.

  • keep_composition (bool, defulat False) – Whether the elemental composition should be the same as in the parents.

get_new_individual(parents)[source]

Function that returns a new individual. Overwrite in subclass.

Comparators meant to be used in symmetry-constrained genetic algorithm (SCGA) for alloy particles.

class acat.ga.symmetric_particle_comparators.ShellSizeComparator(shells, elements=None)[source]

Bases: object

For each given element, compares the sorted sizes of the shells that have the given element. The shells are defined by the symmetry of the particle. Returns True if the sizes are the same, False otherwise. Self-symmetry is considered.

Parameters
  • shells (list of lists) – The atom indices in each shell divided by symmetry. Can be obtained by acat.build.orderings.SymmetricOrderingGenerator.

  • elements (list of strs, default None) – The metal elements of the nanoalloy. Only take into account the elements specified in this list. Default is to take all elements into account.

looks_like(a1, a2)[source]

Return if structure a1 or a2 are similar or not.

class acat.ga.symmetric_particle_comparators.ShellCompositionComparator(shells, elements=None, tol=0)[source]

Bases: object

Compares the elemental compositions of all shells defined by the symmetry of the particle. Returns True if the numbers are the same, False otherwise. Self-symmetry is not considered.

Parameters
  • shells (list of lists) – The atom indices in each shell divided by symmetry. Can be obtained by acat.build.orderings.SymmetricOrderingGenerator.

  • elements (list of strs, default None) – The metal elements of the nanoalloy. Only take into account the elements specified in this list. Default is to take all elements into account.

  • tol (int, default 0) – The maximum number of shells with different elements that two structures are still considered to be look alike.

looks_like(a1, a2)[source]

Return if structure a1 or a2 are similar or not.

Example

All the symmetric particle operators can be easily used with other ASE operators.

As an example we will find the convex hull of ternary NixPtyAu405-x-y truncated octahedral nanoalloys using the ASAP EMT calculator.

The script for a parallel symmetry-constrained genetic algorithm (SCGA) looks as follows:

from acat.build.ordering import SymmetricOrderingGenerator as SOG
from acat.ga.symmetric_particle_operators import SymmetricSubstitute
from acat.ga.symmetric_particle_operators import SymmetricPermutation
from acat.ga.symmetric_particle_operators import SymmetricCrossover
from acat.ga.symmetric_particle_comparators import ShellSizeComparator
from acat.ga.symmetric_particle_comparators import ShellCompositionComparator
from ase.ga.particle_comparator import NNMatComparator
from ase.ga.standard_comparators import SequentialComparator
from ase.ga.offspring_creator import OperationSelector
from ase.ga.population import Population, RankFitnessPopulation
from ase.ga.convergence import GenerationRepetitionConvergence
from ase.ga.data import DataConnection, PrepareDB
from ase.io import read, write
from ase.cluster import Octahedron
from ase.optimize import BFGS
from asap3 import EMT as asapEMT
from multiprocessing import Pool
import os

# Define population.
# Recommand to choose a number that is a multiple of the number of cpu
pop_size = 100

# Generate 100 truncated ocatahedral NixPtyAu405-x-y nanoalloys with
# reflective-circular symmetry. Get the shells at the same time.
particle = Octahedron('Ni', length=9, cutoff=3)
particle.center(vacuum=5.)
sog = SOG(particle, elements=['Ni', 'Pt', 'Au'],
          symmetry='reflective_circular',
          trajectory='starting_generation.traj')
sog.run(max_gen=pop_size, mode='stochastic', verbose=True)
shells = sog.get_shells()
images = read('starting_generation.traj', index=':')

# Instantiate the db
db_name = 'ridge_reflective_circular_NiPtAu_TO405.db'

db = PrepareDB(db_name, cell=particle.cell, population_size=pop_size)

for atoms in images:
    if 'data' not in atoms.info:
        atoms.info['data'] = {}
    db.add_unrelaxed_candidate(atoms, data=atoms.info['data'])

# Connect to the db
db = DataConnection(db_name)

# Define operators
soclist = ([2, 2, 3],
           [SymmetricSubstitute(shells, elements=['Ni', 'Pt', 'Au'], num_muts=3),
            SymmetricPermutation(shells, elements=['Ni', 'Pt', 'Au'], num_muts=1),
            SymmetricCrossover(shells, elements=['Ni', 'Pt', 'Au']),])

op_selector = OperationSelector(*soclist)

# Define comparators
comp = SequentialComparator([ShellSizeComparator(shells, ['Ni', 'Pt', 'Au']),
                             NNMatComparator(0.2, ['Ni', 'Pt', 'Au'])],
                            [0.5, 0.5])

def vf(atoms):
    """Returns the descriptor that distinguishes candidates in the
    niched population."""
    return atoms.get_chemical_formula(mode='hill')

# Give fittest candidates at different compositions equal fitness.
# Use this to find global minima at each composition
pop = RankFitnessPopulation(data_connection=db,
                            population_size=pop_size,
                            comparator=comp,
                            variable_function=vf,
                            exp_function=True,
                            logfile='log.txt')

# Normal fitness ranking irrespective of coverage.
# Use this to find global minimum irrespective of composition
#pop = Population(data_connection=db,
#                 population_size=pop_size,
#                 comparator=comp,
#                 logfile='log.txt')

# Set convergence criteria
cc = GenerationRepetitionConvergence(pop, 5)

# Calculate the relaxed energies for pure Ni405, Pt405 and Au405
pure_pots = {'Ni': 147.532, 'Pt':  86.892, 'Au': 63.566}

# Define the relax function
def relax(atoms, single_point=False):
    atoms.center(vacuum=5.)
    atoms.calc = asapEMT()
    if not single_point:
        opt = BFGS(atoms, logfile=None)
        opt.run(fmax=0.1)
    Epot = atoms.get_potential_energy()
    atoms.info['key_value_pairs']['potential_energy'] = Epot

    # There is a known issue of asapEMT in GA. You can either detach
    # the calculator or re-assign to a SinglePointCalculator
    atoms.calc = None

    # Calculate mixing energy
    syms = atoms.get_chemical_symbols()
    for a in set(syms):
        Epot -= (pure_pots[a] / len(atoms)) * syms.count(a)
    atoms.info['key_value_pairs']['raw_score'] = -Epot

    return atoms

# Relax starting generation
def relax_an_unrelaxed_candidate(atoms):
    if 'data' not in atoms.info:
        atoms.info['data'] = {}
    nncomp = atoms.get_chemical_formula(mode='hill')
    print('Relaxing ' + nncomp)
    relax(atoms)
    db.add_relaxed_step(atoms)

# Create a multiprocessing Pool
pool = Pool(os.cpu_count())
# Perform relaxations in parallel. Especially
# useful when running GA on large nanoparticles
pool.map(relax_an_unrelaxed_candidate, db.get_all_unrelaxed_candidates())
pop.update()

# Number of generations
num_gens = 1000

# Below is the iterative part of the algorithm
gen_num = db.get_generation_number()
for i in range(num_gens):
    # Check if converged
    if cc.converged():
        print('Converged')
        break
    print('Creating and evaluating generation {0}'.format(gen_num + i))

    # Performing procreations in parallel
    def procreation(x):
        # Select an operator and use it
        op = op_selector.get_operator()
        # Select parents for a new candidate
        p1, p2 = pop.get_two_candidates()
        # Pure and binary candidates are not considered
        if len(set(p1.numbers)) < 3:
            return
        parents = [p1, p2]
        offspring, desc = op.get_new_individual(parents)
        # An operator could return None if an offspring cannot be formed
        # by the chosen parents
        if offspring is None:
            return
        nncomp = offspring.get_chemical_formula(mode='hill')
        print('Relaxing ' + nncomp)
        if 'data' not in offspring.info:
            offspring.info['data'] = {}
        relax(offspring)
        db.add_relaxed_candidate(offspring)

    # Create a multiprocessing Pool
    pool = Pool(os.cpu_count())
    # Perform procreations in parallel. Especially
    # useful when running GA on large nanoparticles
    pool.map(procreation, range(pop_size))

    # update the population to allow new candidates to enter
    pop.update()