{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Surface diffusion energy barriers using the Nudged Elastic Band (NEB) method\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n\nfrom ase.build import add_adsorbate, fcc100\nfrom ase.calculators.emt import EMT\nfrom ase.constraints import FixAtoms\nfrom ase.io import read\nfrom ase.mep import NEB\nfrom ase.optimize import BFGS, QuasiNewton\nfrom ase.parallel import world\nfrom ase.visualize.plot import plot_atoms" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First, set up the initial and final states:\n2x2-Al(001) surface with 3 layers and an\nAu atom adsorbed in a hollow site:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "slab = fcc100('Al', size=(2, 2, 3))\nadd_adsorbate(slab, 'Au', 1.7, 'hollow')\nslab.center(axis=2, vacuum=4.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Make sure the structure is correct:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig, ax = plt.subplots()\nplot_atoms(slab, ax)\nax.set_axis_off()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Fix second and third layers:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "mask = [atom.tag > 1 for atom in slab]\nprint(mask)\nslab.set_constraint(FixAtoms(mask=mask))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Use EMT potential:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "slab.calc = EMT()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Initial state:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "qn = QuasiNewton(slab, trajectory='initial.traj')\nqn.run(fmax=0.05)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Final state:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "slab[-1].x += slab.get_cell()[0, 0] / 2\nqn = QuasiNewton(slab, trajectory='final.traj')\nqn.run(fmax=0.05)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Note

Notice how the tags are used to select the constrained atoms

\n\nNow, do the NEB calculation:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "initial = read('initial.traj')\nfinal = read('final.traj')\n\nconstraint = FixAtoms(mask=[atom.tag > 1 for atom in initial])\n\nimages = [initial]\nfor i in range(3):\n image = initial.copy()\n image.calc = EMT()\n image.set_constraint(constraint)\n images.append(image)\n\nimages.append(final)\n\nneb = NEB(images, method='improvedtangent')\nneb.interpolate()\nqn = BFGS(neb, trajectory='neb.traj')\nqn.run(fmax=0.05)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Visualize the results with::\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig, ax = plt.subplots()\nplot_atoms(final, ax)\nax.set_axis_off()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "or from the command line:\n\n $ ase gui neb.traj@-5:\n\nand select Tools->NEB.\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can also create a series of plots like above, that show the progression\nof the NEB relaxation, directly at the command line::\n\n $ ase nebplot --share-x --share-y neb.traj\n\nFor more customizable analysis of the output of many NEB jobs, you can use\nthe :class:`ase.mep.NEBTools` class. Some examples of its use are below; the\nfinal example was used to make the figure you see above.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n\nfrom ase.io import read\nfrom ase.mep import NEBTools\n\nimages = read('neb.traj@-5:')\n\nnebtools = NEBTools(images)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get the calculated barrier and the energy change of the reaction.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "Ef, dE = nebtools.get_barrier()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get the barrier without any interpolation between highest images.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "Ef, dE = nebtools.get_barrier(fit=False)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get the actual maximum force at this point in the simulation.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "max_force = nebtools.get_fmax()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create a figure like that coming from ASE-GUI.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig = nebtools.plot_band()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create a figure with custom parameters.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig = plt.figure(figsize=(5.5, 4.0))\nax = fig.add_axes((0.15, 0.15, 0.8, 0.75))\nnebtools.plot_band(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Note

For this reaction, the reaction coordinate is very simple: The\n *x*-coordinate of the Au atom. In such cases, the NEB method is\n overkill, and a simple constraint method should be used like in this\n tutorial: `constraints diffusion tutorial`.

\n\n.. seealso::\n\n * :mod:`ase.mep.neb`\n * :mod:`ase.constraints`\n * `constraints diffusion tutorial`\n * :func:`~ase.build.fcc100`\n\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Restarting NEB\n\nRestart NEB from the trajectory file:\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "read the last structures (of 5 images used in NEB)\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "images = read('neb.traj@-5:')\n\nfor i in range(1, len(images) - 1):\n images[i].calc = EMT()\n\nneb = NEB(images, method='improvedtangent')\nqn = BFGS(neb, trajectory='neb_restart.traj')\nqn.run(fmax=0.005)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Parallelizing over images with MPI\n\nInstead of having one process do the calculations for all three\ninternal images in turn, it will be faster to have three processes do\none image each. In order to be able to run python with MPI\nyou need a special parallel python interpreter, for example gpaw python\n(see [GPAW parallel runs](https://gpaw.readthedocs.io/documentation/parallel_runs/parallel_runs.html))\nand set ``parallel=True`` in the NEB calculation.\n\nThe example below can then be run\nwith ``gpaw -P3 python diffusion_parallel.py``:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "initial = read('initial.traj')\nfinal = read('final.traj')\n\nconstraint = FixAtoms(mask=[atom.tag > 1 for atom in initial])\n\nn_images = 1 # Set to number of processes you use with mpiexec\nimages = [initial]\nj = world.rank * n_images // world.size # my image number\nfor i in range(n_images):\n image = initial.copy()\n if i == j:\n image.calc = EMT()\n image.set_constraint(constraint)\n images.append(image)\nimages.append(final)\n\nneb = NEB(images, parallel=True, method='improvedtangent')\nneb.interpolate()\nqn = BFGS(neb, trajectory='neb.traj')\nqn.run(fmax=0.05)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.7" } }, "nbformat": 4, "nbformat_minor": 0 }