
.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples_generated/tutorials/constraint_diffusion.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_examples_generated_tutorials_constraint_diffusion.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_examples_generated_tutorials_constraint_diffusion.py:

.. _constraints diffusion tutorial:

============================================================
Constrained Calculations - Surface diffusion energy barriers
============================================================

In this tutorial, we will calculate the energy barrier that was found
using the :mod:`NEB <ase.mep.neb>` method in the :ref:`diffusion tutorial`
tutorial.  Here, we use a simple :class:`~ase.constraints.FixedPlane`
constraint that forces the Au atom to relax in the *yz*-plane only:

.. GENERATED FROM PYTHON SOURCE LINES 12-24

.. code-block:: Python


    from ase.build import add_adsorbate, fcc100
    from ase.calculators.emt import EMT
    from ase.constraints import FixAtoms, FixedPlane
    from ase.optimize import QuasiNewton

    # 2x2-Al(001) surface with 3 layers and an
    # Au atom adsorbed in a hollow site:
    slab = fcc100('Al', size=(2, 2, 3))
    add_adsorbate(slab, 'Au', 1.7, 'hollow')
    slab.center(axis=2, vacuum=4.0)








.. GENERATED FROM PYTHON SOURCE LINES 25-26

We can visualize the structure with ase visualize:

.. GENERATED FROM PYTHON SOURCE LINES 26-38

.. code-block:: Python

    import matplotlib.pyplot as plt

    from ase.visualize.plot import plot_atoms

    fig, (ax1, ax2) = plt.subplots(1, 2)
    plot_atoms(slab, ax1)
    plot_atoms(slab, ax2, rotation='-90x')
    ax1.set_title('top view')
    ax2.set_title('side view')
    ax1.set_axis_off()
    ax2.set_axis_off()




.. image-sg:: /examples_generated/tutorials/images/sphx_glr_constraint_diffusion_001.png
   :alt: top view, side view
   :srcset: /examples_generated/tutorials/images/sphx_glr_constraint_diffusion_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 39-40

Alternatively, you can use also use view directly:

.. GENERATED FROM PYTHON SOURCE LINES 40-44

.. code-block:: Python


    #   $ from ase.visualize import view
    #   $ view(slab)








.. GENERATED FROM PYTHON SOURCE LINES 45-46

We can now continue fixing the atoms in the slab:

.. GENERATED FROM PYTHON SOURCE LINES 46-57

.. code-block:: Python


    # Fix second and third layers:
    mask = [atom.tag > 1 for atom in slab]
    # print(mask)
    fixlayers = FixAtoms(mask=mask)

    # Constrain the last atom (Au atom) to move only in the yz-plane:
    plane = FixedPlane(-1, (1, 0, 0))

    slab.set_constraint([fixlayers, plane])








.. GENERATED FROM PYTHON SOURCE LINES 58-61

Now we can perform the calculation optimizing the displacement
of the gold atom along the x-axis.
We do structure optimization here using the EMT potential:

.. GENERATED FROM PYTHON SOURCE LINES 61-83

.. code-block:: Python


    # Use EMT potential:
    slab.calc = EMT()

    for i in range(5):
        qn = QuasiNewton(slab, trajectory=f'mep{i}.traj')
        qn.run(fmax=0.05)
        # Move gold atom along x-axis:
        slab[-1].x += slab.get_cell()[0, 0] / 8

    # Let's visualize the saved trajectory.
    # Here is code to visualize
    # a side-view of the path (unit cell repeated twice):


    from ase.io import read

    configs = [read(f'mep{i}.traj', '-1') for i in range(5)]

    # for easier visualization, let's repeat the structures
    configs_repeated = [config.repeat((2, 1, 1)) for config in configs]





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

                    Step[ FC]     Time          Energy          fmax
    BFGSLineSearch:    0[  0] 09:17:31        3.323870       0.2462
    BFGSLineSearch:    1[  1] 09:17:31        3.314754       0.0378
                    Step[ FC]     Time          Energy          fmax
    BFGSLineSearch:    0[  0] 09:17:31        3.730166       1.6460
    BFGSLineSearch:    1[  2] 09:17:31        3.546230       0.2605
    BFGSLineSearch:    2[  4] 09:17:31        3.513795       0.1080
    BFGSLineSearch:    3[  6] 09:17:31        3.502279       0.0500
                    Step[ FC]     Time          Energy          fmax
    BFGSLineSearch:    0[  0] 09:17:31        3.764262       1.0444
    BFGSLineSearch:    1[  2] 09:17:31        3.708446       0.1203
    BFGSLineSearch:    2[  4] 09:17:31        3.695835       0.0668
    BFGSLineSearch:    3[  6] 09:17:31        3.691660       0.0618
    BFGSLineSearch:    4[  7] 09:17:31        3.690089       0.0300
                    Step[ FC]     Time          Energy          fmax
    BFGSLineSearch:    0[  0] 09:17:31        3.580034       0.6521
    BFGSLineSearch:    1[  1] 09:17:31        3.544894       0.1627
    BFGSLineSearch:    2[  3] 09:17:31        3.516850       0.1946
    BFGSLineSearch:    3[  4] 09:17:31        3.508545       0.1325
    BFGSLineSearch:    4[  5] 09:17:31        3.503316       0.1022
    BFGSLineSearch:    5[  6] 09:17:31        3.500392       0.0550
    BFGSLineSearch:    6[  8] 09:17:31        3.499584       0.0392
                    Step[ FC]     Time          Energy          fmax
    BFGSLineSearch:    0[  0] 09:17:31        3.463634       0.5332
    BFGSLineSearch:    1[  2] 09:17:31        3.356702       0.4028
    BFGSLineSearch:    2[  3] 09:17:31        3.333298       0.1536
    BFGSLineSearch:    3[  4] 09:17:32        3.322604       0.1003
    BFGSLineSearch:    4[  5] 09:17:32        3.318477       0.0658
    BFGSLineSearch:    5[  6] 09:17:32        3.316534       0.0483




.. GENERATED FROM PYTHON SOURCE LINES 84-85

We can visualize the structures with `ase.visualize.plot.animate`:

.. GENERATED FROM PYTHON SOURCE LINES 85-95

.. code-block:: Python


    from ase.visualize.plot import animate

    animate(
        configs_repeated,
        ax=None,
        interval=500,  # in ms; same default value as in FuncAnimation
        rotation=('-90x,0y,0z'),
    )




.. video:: /examples_generated/tutorials/images/sphx_glr_constraint_diffusion_002.mp4
   :class: sphx-glr-single-img
   :height: 480
   :width: 640
   :autoplay:
   :loop:



.. rst-class:: sphx-glr-script-out

 .. code-block:: none


    <matplotlib.animation.FuncAnimation object at 0x7ced210b4ec0>



.. GENERATED FROM PYTHON SOURCE LINES 96-97

Let's plot the energy and look at the barrier.

.. GENERATED FROM PYTHON SOURCE LINES 97-107

.. code-block:: Python


    # get the potential energies of the structures
    energies = [config.get_potential_energy() for config in configs]
    # set last energy value to 0 for easier comparison
    energies = [energy - energies[-1] for energy in energies]
    plt.ylabel(r'$E(i) - E_{\mathrm{final}}$ (meV)')
    plt.xlabel('Image number')
    plt.plot(range(1, len(energies) + 1), energies)





.. image-sg:: /examples_generated/tutorials/images/sphx_glr_constraint_diffusion_003.png
   :alt: constraint diffusion
   :srcset: /examples_generated/tutorials/images/sphx_glr_constraint_diffusion_003.png
   :class: sphx-glr-single-img


.. rst-class:: sphx-glr-script-out

 .. code-block:: none


    [<matplotlib.lines.Line2D object at 0x7ced210b67b0>]



.. GENERATED FROM PYTHON SOURCE LINES 108-115

The barrier is found to
be 0.35 eV - exactly as in the :ref:`NEB <diffusion tutorial>`
tutorial.

The result can also be analysed with the
command :command:`ase gui mep?.traj -n
-1` (choose :menuselection:`Tools --> NEB`).

.. GENERATED FROM PYTHON SOURCE LINES 118-124

.. seealso::

  * :mod:`ase.mep.neb`
  * :mod:`ase.constraints`
  * :ref:`diffusion tutorial`
  * :func:`~ase.build.fcc100`


.. _sphx_glr_download_examples_generated_tutorials_constraint_diffusion.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: constraint_diffusion.ipynb <constraint_diffusion.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: constraint_diffusion.py <constraint_diffusion.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: constraint_diffusion.zip <constraint_diffusion.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_
