Genetic algorithm¶
Optimize adsorbate overlayer 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, species_probabilities=None, 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.
-
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.
-
classmethod
-
class
acat.ga.adsorbate_operators.
AddAdsorbate
(adsorbate_species, species_probabilities=None, 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, adsorption_sites=None, surface=None, allow_6fold=False, site_preference=None, surface_preference=None, max_coverage=None, tilt_angle=None, num_muts=1, dmax=3.0)[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.
species_probabilities (dict, default None) – A dictionary that contains keys of each adsorbate species and values of their probabilities of adding onto the surface. Adding adsorbate species with equal probability if not specified.
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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
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.
site_preference (str of list of strs, defualt None) – The site type(s) 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.
max_coverage (float, default None) – The maximum allowed coverage on the surface. Coverage is defined as (number of surface occupied sites / number of surface atoms). The maximum coverage is solely governed by min_adsorbate_distance if max_coverage is not specified.
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 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
class
acat.ga.adsorbate_operators.
RemoveAdsorbate
(adsorbate_species, adsorption_sites=None, surface=None, allow_6fold=False, site_preference=None, surface_preference=None, min_coverage=None, num_muts=1, dmax=3.0)[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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
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.
site_preference (str or list of strs, defualt None) – The site type(s) that has higher priority to detach adsorbates.
surface_preference (str, default None) – The surface type that has higher priority to detach adsorbates. Please only use this for nanoparticles.
min_coverage (float, default None) – The minimum allowed coverage on the surface. Coverage is defined as (number of surface occupied sites / number of surface atoms). The minimum coverage is 0 if min_coverage is not specified.
num_muts (int, default 1) – The number of times to perform this operation.
dmax (float, default 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
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, adsorption_sites=None, surface=None, allow_6fold=False, site_preference_from=None, surface_preference_from=None, site_preference_to=None, surface_preference_to=None, tilt_angle=None, num_muts=1, dmax=3.0)[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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
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.
site_preference_from (str or list of strs, defualt None) – The site type(s) 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 or list of strs, defualt None) – The site type(s) 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 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
class
acat.ga.adsorbate_operators.
ReplaceAdsorbate
(adsorbate_species, species_probabilities=None, 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, adsorption_sites=None, surface=None, allow_6fold=False, site_preference=None, surface_preference=None, tilt_angle=None, num_muts=1, dmax=3.0)[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.
species_probabilities (dict, default None) – A dictionary that contains keys of each adsorbate species and values of their probabilities of replacing an adsorbate on the surface. Adding adsorbate species with equal probability if not specified.
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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
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.
site_preference (str or list of strs, defualt None) – The site type(s) that has higher priority to replace adsorbates.
surface_preference (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 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
class
acat.ga.adsorbate_operators.
CutSpliceCrossoverWithAdsorbates
(adsorbate_species, blmin, 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}, keep_composition=True, fix_coverage=False, min_adsorbate_distance=2.0, allow_6fold=False, rotate_vectors=None, rotate_angles=None, dmax=3.0)[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 (only returns one of them). The indexing of the atoms is not preserved. Please only use this operator if the particle is allowed to change shape.
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}
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.
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.
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 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
class
acat.ga.adsorbate_operators.
SimpleCutSpliceCrossoverWithAdsorbates
(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}, keep_composition=True, fix_coverage=False, min_adsorbate_distance=2.0, adsorption_sites=None, allow_6fold=False, dmax=3.0)[source]¶ Bases:
acat.ga.adsorbate_operators.AdsorbateOperator
Crossover that divides two particles through a plane in space and merges the symbols of two halves from different particles with adosrbates together (only returns one of them). The indexing of the atoms is preserved. Please only use this operator with other operators that also preserves the indexing.
It keeps the correct composition by randomly assigning elements in the new particle.
- 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.
keep_composition (bool, default True) – Boolean that signifies if the composition should 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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
allow_6fold (bool, default False) – Whether to allow the adsorption on 6-fold subsurf sites underneath fcc hollow sites.
dmax (float, default 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
-
class
acat.ga.adsorbate_operators.
AdsorbateCatalystCrossover
(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}, adsorption_sites=None, surface=None, allow_6fold=False, dmax=3.0)[source]¶ Bases:
acat.ga.adsorbate_operators.AdsorbateOperator
Crossover that divides two particles or two slabs by the adsorbate -catalyst interfaces and exchange all adsorbates (only returns one of them). The indexing of the atoms is preserved. Please only use this operator with other operators that also preserves the indexing.
The composition or the coverage is fixed if it is preserved by all other operators being used.
- 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.
adsorption_sites (acat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object, default None) – Provide the built-in adsorption sites class to accelerate the genetic algorithm. Make sure all the operators used with this operator preserve the indexing of the atoms.
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.
dmax (float, default 3.) – The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
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
-
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)
-
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:
-
class
acat.ga.adsorbate_comparators.
AdsorptionGraphComparator
(adsorption_sites, composition_effect=True, fragmentation=True, subsurf_effect=False, full_effect=False, dmax=3.0)[source]¶ Bases:
object
Compares the graph of adsorbate overlayer + surface atoms and returns True if they are isomorphic with node matches. Before checking graph isomorphism, a cheap label match is used to reject graphs that are impossible to be isomorphic.
The graphs can be quite costly to obtain every time a graph is required (and disk intensive if saved), thus it makes sense to get the graph along with e.g. the potential energy and save it in atoms.info[‘data’][‘graph’].
Parameters:
- adsorption_sitesacat.adsorption_sites.ClusterAdsorptionSites object or acat.adsorption_sites.SlabAdsorptionSites object
Provide the acat built-in adsorption sites class to accelerate the pattern generation. Make sure all the structures have the same atom indexing.
- composition_effectbool, default True
Whether to consider sites with different elemental compositions as different sites. It is recommended to set composition_effet=False for monometallics.
- subsurf_effectbool, default False
Whether to take subsurface atoms into consideration when checking uniqueness. Could be important for surfaces like fcc100.
- full_effectbool, default False
Take the whole catalyst into consideration when generating graph.
- dmaxfloat, default 3.
The maximum bond length (in Angstrom) between the site and the bonding atom that should be considered as an adsorbate.
- fragmentationbool, default True
Whether to cut multidentate species into fragments. This ensures that multidentate species with different orientations are considered as different adlayer patterns.
Example1
All the adsorbate operators and comparators can be easily used with other operators and comparators. AddAdsorbate
, RemoveAdsorbate
, MoveAdsorbate
, ReplaceAdsorbate
and AdsorbateCatalystCrossover
operators can be used for both non-periodic nanoparticles and periodic surface slabs. CutSpliceCrossoverWithAdsorbates
and SimpleCutSpliceCrossoverWithAdsorbates
operators only work for nanoparticles, and the latter is recommended. To accelerate the GA, provide adsorsption sites and use indexing-preserved operators implemented in ACAT.
As an example we will simultaneously optimize both the adsorbate overlayer 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.adsorption_sites import ClusterAdsorptionSites
from acat.adsorbate_coverage import ClusterAdsorbateCoverage
from acat.build.ordering import RandomOrderingGenerator as ROG
from acat.build.adlayer import min_dist_coverage_pattern
from acat.ga.adsorbate_operators import (AddAdsorbate, RemoveAdsorbate,
MoveAdsorbate, ReplaceAdsorbate,
SimpleCutSpliceCrossoverWithAdsorbates)
# Import particle_mutations from acat instead of ase to get the indexing-preserved version
from acat.ga.particle_mutations import (RandomPermutation, COM2surfPermutation,
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, get_nnmat
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 numpy as np
import time
import os
# Define population
# Recommend 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(num_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 = min_dist_coverage_pattern(atoms, adsorbate_species=species,
min_adsorbate_distance=dmin)
patterns.append(pattern)
# Get the adsorption sites. Composition does not matter in GA
sas = ClusterAdsorptionSites(particle, composition_effect=False)
# 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, 2],
[Rich2poorPermutation(elements=['Ni', 'Pt'], num_muts=5),
Poor2richPermutation(elements=['Ni', 'Pt'], num_muts=5),
RandomPermutation(elements=['Ni', 'Pt'], num_muts=5),
AddAdsorbate(species, adsorption_sites=sas, num_muts=5),
RemoveAdsorbate(species, adsorption_sites=sas, num_muts=5),
MoveAdsorbate(species, adsorption_sites=sas, num_muts=5),
ReplaceAdsorbate(species, adsorption_sites=sas, num_muts=5),
SimpleCutSpliceCrossoverWithAdsorbates(species, keep_composition=True,
adsorption_sites=sas),])
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)
time.sleep(1) # Add delay to avoid stack overflow (only for testing)
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
# Parallelize nnmat calculations to accelerate NNMatComparator
atoms.info['data']['nnmat'] = get_nnmat(atoms)
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 only for testing
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())
pool.close()
pool.join()
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()
# Assign rng with a random seed
np.random.seed(random.randint(1, 10000))
pop.rng = np.random
# 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:
return
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, single_point=True) # Single point only for testing
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))
pool.close()
pool.join()
# update the population to allow new candidates to enter
pop.update()
Example2
If we want to study the same system but fix the adsorbate coverage to be exactly 20 adsorbates, we should only use MoveAdsorbate
, ReplaceAdsorbate
and SimpleCutSpliceCrossoverWithAdsorbates
operators with same particle operators. Remember to set fix_coverage to True in SimpleCutSpliceCrossoverWithAdsorbates
.
The script for a fixed-coverage parallel genetic algorithm now looks as follows:
from acat.settings import adsorbate_elements
from acat.adsorption_sites import ClusterAdsorptionSites
from acat.adsorbate_coverage import ClusterAdsorbateCoverage
from acat.build.ordering import RandomOrderingGenerator as ROG
from acat.build.adlayer import StochasticPatternGenerator as SPG
from acat.ga.adsorbate_operators import (AddAdsorbate, RemoveAdsorbate,
MoveAdsorbate, ReplaceAdsorbate,
SimpleCutSpliceCrossoverWithAdsorbates)
# Import particle_mutations from acat instead of ase to get the indexing-preserved version
from acat.ga.particle_mutations import (RandomPermutation, COM2surfPermutation,
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, get_nnmat
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 numpy as np
import time
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(num_gen=pop_size)
# Generate random coverage patterns of 20 adsorbates on nanoparticles
species=['H', 'C', 'O', 'OH', 'CO', 'CH', 'CH2', 'CH3']
images = read('starting_generation.traj', index=':')
num_ads = 20 # Number of adsorbates, fix at this coverage throughout GA
for _ in range(num_ads): # This will take quite some time
spg = SPG(images, adsorbate_species=['CO','OH','N'],
min_adsorbate_distance=3.,
composition_effect=True)
spg.run(num_gen=pop_size, action='add', unique=False)
images = read('patterns.traj')
patterns = read('patterns.traj', index=':')
# Get the adsorption sites. Composition does not matter in GA
sas = ClusterAdsorptionSites(particle, composition_effect=False)
# 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, 2],
[Rich2poorPermutation(elements=['Ni', 'Pt'], num_muts=5),
Poor2richPermutation(elements=['Ni', 'Pt'], num_muts=5),
RandomPermutation(elements=['Ni', 'Pt'], num_muts=5),
MoveAdsorbate(species, adsorption_sites=sas, num_muts=5),
ReplaceAdsorbate(species, adsorption_sites=sas, num_muts=5),
SimpleCutSpliceCrossoverWithAdsorbates(species, keep_composition=True,
fix_coverage=True,
adsorption_sites=sas),])
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)
time.sleep(1) # Add delay to avoid stack overflow (just for testing)
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
# Parallelize nnmat calculations to accelerate NNMatComparator
atoms.info['data']['nnmat'] = get_nnmat(atoms)
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 only for testing
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())
pool.close()
pool.join()
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()
# Assign rng with a random seed
np.random.seed(random.randint(1, 10000))
pop.rng = np.random
# 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:
return
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, single_point=True) # Single point only for testing
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))
pool.close()
pool.join()
# 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).
-
class
acat.ga.group_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.group_operators.
GroupSubstitute
(groups, elements=None, num_muts=1)[source]¶ Bases:
acat.ga.group_operators.Mutation
Substitute all the atoms in a group with a different element. The elemental composition cannot be fixed.
- Parameters
groups (list of lists) – The atom indices in each user-divided group. Can be obtained by acat.build.ordering.SymmetricClusterOrderingGenerator or acat.build.ordering.OrderedSlabOrderingGenerator.
elements (list of strs, default None) – 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.
-
class
acat.ga.group_operators.
GroupPermutation
(groups, elements=None, keep_composition=False, num_muts=1)[source]¶ Bases:
acat.ga.group_operators.Mutation
Permutes the elements in two random groups. The elemental composition can be fixed.
- Parameters
groups (list of lists) – The atom indices in each user-divided group. Can be obtained by acat.build.ordering.SymmetricClusterOrderingGenerator or acat.build.ordering.OrderedSlabOrderingGenerator.
elements (list of strs, default None) – 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.
-
class
acat.ga.group_operators.
Crossover
[source]¶ Bases:
ase.ga.offspring_creator.OffspringCreator
Base class for all particle crossovers. Do not call this class directly.
-
class
acat.ga.group_operators.
GroupCrossover
(groups, elements=None, keep_composition=False)[source]¶ Bases:
acat.ga.group_operators.Crossover
Merge the elemental distributions in two half groups from different structures together. The elemental composition can be fixed.
- Parameters
groups (list of lists) – The atom indices in each user-divided group. Can be obtained by acat.build.ordering.SymmetricClusterOrderingGenerator or acat.build.ordering.OrderedSlabOrderingGenerator.
elements (list of strs, default None) – 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.
Comparators meant to be used in symmetry-constrained genetic algorithm (SCGA).
-
class
acat.ga.group_comparators.
GroupSizeComparator
(groups, elements=None)[source]¶ Bases:
object
For each given element, compares the sorted sizes of the user-divided groups that have the given element. Returns True if the sizes are the same, False otherwise. Self-symmetry is considered for particles.
- Parameters
groups (list of lists) – The atom indices in each user-divided group. Can be obtained by acat.build.ordering.SymmetricClusterOrderingGenerator or acat.build.ordering.OrderedSlabOrderingGenerator.
elements (list of strs, default None) – Only take into account the elements specified in this list. Default is to take all elements into account.
-
class
acat.ga.group_comparators.
GroupCompositionComparator
(groups, elements=None, tol=0)[source]¶ Bases:
object
Compares the elemental compositions of all user-divided groups. Returns True if the numbers are the same, False otherwise. Self-symmetry is not considered for particles.
- Parameters
groups (list of lists) – The atom indices in each user-divided group. Can be obtained by acat.build.ordering.SymmetricClusterOrderingGenerator or acat.build.ordering.OrderedSlabOrderingGenerator.
elements (list of strs, default None) – 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 groups with different elements that two structures are still considered to be look alike.
Example1
All the group operators and comparators can be easily used with other indexing-preserved operators and comparators. All operators can be used for both non-periodic nanoparticles and periodic surface slabs.
As an example we will find the convex hull of ternary NixPtyAu405-x-y truncated octahedral nanoalloys using the ASAP EMT calculator. Note that we must first align the symmetry axis of interest to the z direction. Here we want to study the 3-fold mirror circular symmetry around the C3 axis of the particle.
The script for a parallel symmetry-constrained genetic algorithm (SCGA) looks as follows:
from acat.build.ordering import SymmetricClusterOrderingGenerator as SCOG
from acat.ga.group_operators import (GroupSubstitute,
GroupPermutation,
GroupCrossover)
from acat.ga.group_comparators import (GroupSizeComparator,
GroupCompositionComparator)
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.ga.utilities import get_nnmat
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 numpy as np
import os
# Define population.
# Recommend to choose a number that is a multiple of the number of cpu
pop_size = 100
# Build and rotate the particle so that C3 axis is aligned to z direction
particle = Octahedron('Ni', length=9, cutoff=3)
particle.rotate(45, 'x')
particle.rotate(-35.29, 'y')
particle.center(vacuum=5.)
# Generate 100 truncated ocatahedral NixPtyAu405-x-y nanoalloys with
# mirror circular symmetry. Get the groups at the same time.
scog = SCOG(particle, elements=['Ni', 'Pt', 'Au'],
symmetry='mirror_circular',
trajectory='starting_generation.traj')
scog.run(max_gen=pop_size, mode='stochastic', verbose=True)
groups = scog.get_groups()
images = read('starting_generation.traj', index=':')
# Instantiate the db
db_name = 'ridge_mirror_circular_NiPtAu_TO405_C3.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],
[GroupSubstitute(groups, elements=['Ni', 'Pt', 'Au'], num_muts=3),
GroupPermutation(groups, elements=['Ni', 'Pt', 'Au'], num_muts=1),
GroupCrossover(groups, elements=['Ni', 'Pt', 'Au']),])
op_selector = OperationSelector(*soclist)
# Define comparators
comp = SequentialComparator([GroupSizeComparator(groups, ['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
# Parallelize nnmat calculations to accelerate NNMatComparator
atoms.info['data']['nnmat'] = get_nnmat(atoms)
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())
pool.close()
pool.join()
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()
# Assign rng with a random seed
np.random.seed(random.randint(1, 10000))
pop.rng = np.random
# 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))
pool.close()
pool.join()
# update the population to allow new candidates to enter
pop.update()
Example2
If we want to study the same system but target 3 compositions: Ni0.5Pt0.25Au0.25, Ni0.25Pt0.5Au0.25 and Ni0.25Pt0.25Au0.5, we should not use GroupSubstitute
operator and set keep_composition to True in GroupPermutation
and GroupCrossover
operators. The tolerance of the intitial compositions can be controlled by the eps parameter in SymmetricClusterOrderingGenerator.run
.
The script for a fixed-composition parallel genetic algorithm now looks as follows:
from acat.build.ordering import SymmetricClusterOrderingGenerator as SCOG
from acat.ga.group_operators import (GroupSubstitute,
GroupPermutation,
GroupCrossover)
from acat.ga.group_comparators import (GroupSizeComparator,
GroupCompositionComparator)
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.ga.utilities import get_nnmat
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 numpy as np
import os
# Define population.
# Recommend to choose a number that is a multiple of the number of cpu
pop_size = 100
# Build and rotate the particle so that C3 axis is aligned to z direction
particle = Octahedron('Ni', length=9, cutoff=3)
particle.rotate(45, 'x')
particle.rotate(-35.29, 'y')
particle.center(vacuum=5.)
# Generate 100 truncated ocatahedral NixPtyAu405-x-y nanoalloys with
# mirror circular symmetry and concentrations of {x=0.5, y=0.25},
# {x=0.25, y=0.5} and {x=y=0.25}. The concentrations are allowed to
# vary by a range of 2*eps=0.1. Get the groups at the same time.
scog1 = SCOG(particle, elements=['Ni', 'Pt'],
symmetry='mirror_circular',
cutoff=1.,
secondary_symmetry='chemical',
secondary_cutoff=.1,
composition={'Ni': 0.5, 'Pt': 0.25, 'Au': 0.25},
trajectory='starting_generation.traj')
groups = scog1.get_groups()
scog1.run(max_gen=33, mode='stochastic', eps=0.05, verbose=True)
scog2 = SCOG(particle, elements=['Ni', 'Pt'],
symmetry='mirror_circular',
cutoff=1.,
secondary_symmetry='chemical',
secondary_cutoff=.1,
composition={'Ni': 0.25, 'Pt': 0.5, 'Au': 0.25},
trajectory='starting_generation.traj',
append_trajectory=True)
scog2.run(max_gen=33, mode='stochastic', eps=0.05, verbose=True)
scog3 = SCOG(particle, elements=['Ni', 'Pt'],
symmetry='mirror_circular',
cutoff=1.,
secondary_symmetry='chemical',
secondary_cutoff=.1,
composition={'Ni': 0.25, 'Pt': 0.25, 'Au': 0.5},
trajectory='starting_generation.traj',
append_trajectory=True)
scog3.run(max_gen=34, mode='stochastic', eps=0.05, verbose=True)
images = read('starting_generation.traj', index=':')
# Instantiate the db
db_name = 'ridge_mirror_circular_NiPtAu_TO405_C3.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, now set keep_composition=True
# GroupSubstitute cannot keep the composition so it's not used
soclist = ([4, 6],
[GroupPermutation(groups, elements=['Ni', 'Pt', 'Au'],
keep_composition=True, num_muts=1),
GroupCrossover(groups, elements=['Ni', 'Pt', 'Au'],
keep_composition=True),])
op_selector = OperationSelector(*soclist)
# Define comparators
comp = SequentialComparator([GroupSizeComparator(groups, ['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
# Parallelize nnmat calculations to accelerate NNMatComparator
atoms.info['data']['nnmat'] = get_nnmat(atoms)
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())
pool.close()
pool.join()
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()
# Assign rng with a random seed
np.random.seed(random.randint(1, 10000))
pop.rng = np.random
# 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))
pool.close()
pool.join()
# update the population to allow new candidates to enter
pop.update()