.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples_generated/tutorials/bulk.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_generated_tutorials_bulk.py: .. _bulk: Bulk Structures and Relaxations =============================== Here, we create bulk structures and optimize them to their ideal bulk properties .. GENERATED FROM PYTHON SOURCE LINES 10-22 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from ase.build import bulk from ase.calculators.emt import EMT from ase.eos import EquationOfState from ase.filters import FrechetCellFilter from ase.optimize import BFGS from ase.units import kJ from ase.visualize.plot import plot_atoms .. GENERATED FROM PYTHON SOURCE LINES 23-44 Setting up bulk structures -------------------------- ASE provides three frameworks for setting up bulk structures: * :func:`ase.build.bulk`. Knows lattice types and lattice constants for elemental bulk structures and a few compounds, but with limited customization. * :func:`ase.spacegroup.crystal`. Creates atoms from typical crystallographic information such as spacegroup, lattice parameters, and basis. * :mod:`ase.lattice`. Creates atoms explicitly from lattice and basis. Let's run a simple bulk calculation. We use :func:`ase.build.bulk` to get a primitive cell of silver, and then visualize it. Silver is known to form an FCC structure, so presumably the function returned a primitive FCC cell. You can, e.g., use the ASE GUI to repeat the structure and recognize the A-B-C stacking. .. GENERATED FROM PYTHON SOURCE LINES 45-57 .. code-block:: Python atoms = bulk('Ag') fig, ax = plt.subplots() plot_atoms(atoms * (3, 3, 3), ax=ax) ax.set_xlabel(r'$x (\AA)$') ax.set_ylabel(r'$y (\AA)$') fig.tight_layout() # For interactive use of the ASE GUI: # view(atoms * (3,3,3)) .. image-sg:: /examples_generated/tutorials/images/sphx_glr_bulk_001.png :alt: bulk :srcset: /examples_generated/tutorials/images/sphx_glr_bulk_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 58-60 ASE should also be able to verify that it really is a primitive FCC cell and tell us what lattice constant was chosen: .. GENERATED FROM PYTHON SOURCE LINES 60-63 .. code-block:: Python print(f'Bravais lattice: {atoms.cell.get_bravais_lattice()}') .. rst-class:: sphx-glr-script-out .. code-block:: none Bravais lattice: FCC(a=4.09) .. GENERATED FROM PYTHON SOURCE LINES 64-73 Periodicity ----------- Periodic structures in ASE are represented using ``atoms.cell`` and ``atoms.pbc``. * The cell is a :class:`~ase.cell.Cell` object which represents the crystal lattice with three vectors. * ``pbc`` is an array of three booleans indicating whether the system is periodic in each direction. .. GENERATED FROM PYTHON SOURCE LINES 74-78 .. code-block:: Python print('Cell:\n', atoms.cell.round(3)) print('Periodicity: ', atoms.pbc) .. rst-class:: sphx-glr-script-out .. code-block:: none Cell: [[0. 2.045 2.045] [2.045 0. 2.045] [2.045 2.045 0. ]] Periodicity: [ True True True] .. GENERATED FROM PYTHON SOURCE LINES 79-93 Equation of state ----------------- We can find the optimal lattice parameter and calculate the bulk modulus by doing an equation-of-state calculation. This means sampling the energy and lattice constant over a range of values to get the minimum as well as the curvature, which gives us the bulk modulus. The online ASE docs already provide a tutorial on how to perform equation-of-state calculations: (:ref:`eos_example`) First, we calculate the volume and potential energy, whilst scaling the atoms' cell. Here, we use ase's empirical calculator ``EMT``: .. GENERATED FROM PYTHON SOURCE LINES 94-108 .. code-block:: Python calc = EMT() cell = atoms.get_cell() volumes = [] energies = [] for x in np.linspace(0.95, 1.05, 5): atoms_copy = atoms.copy() atoms_copy.calc = calc atoms_copy.set_cell(cell * x, scale_atoms=True) atoms_copy.get_potential_energy() volumes.append(atoms_copy.get_volume()) energies.append(atoms_copy.get_potential_energy()) .. GENERATED FROM PYTHON SOURCE LINES 109-112 Then, via :func:`ase.eos.EquationOfState`, we can calculate and plot the bulk modulus: .. GENERATED FROM PYTHON SOURCE LINES 112-125 .. code-block:: Python eos = EquationOfState(volumes, energies) v0, e0, B = eos.fit() ax = eos.plot() ax.axhline(e0, linestyle='--', alpha=0.5, color='black') ax.axvline(v0, linestyle='--', alpha=0.5, color='red') print(f'Minimum Volume = {v0:.3f}AA^3') print(f'Minimum Energy = {e0:.3f}eV') print(f'Bulk modulus = {B / kJ * 1.0e24:.3f} GPa') plt.show() .. image-sg:: /examples_generated/tutorials/images/sphx_glr_bulk_002.png :alt: sj: E: -0.000 eV, V: 16.776 Å$^3$, B: 100.124 GPa :srcset: /examples_generated/tutorials/images/sphx_glr_bulk_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Minimum Volume = 16.776AA^3 Minimum Energy = -0.000eV Bulk modulus = 100.124 GPa .. GENERATED FROM PYTHON SOURCE LINES 126-134 Bulk Optimization ----------------- We can also find the optimal, relaxed cell via variable-cell relaxation. This requires setting a filter, in this case the :func:`ase.filters.FrechetCellFilter` filter, which allows minimizing both the atomic forces and the unit cell stress. .. GENERATED FROM PYTHON SOURCE LINES 135-146 .. code-block:: Python original_lattice = atoms.cell.get_bravais_lattice() calc = EMT() atoms.calc = calc opt = BFGS(FrechetCellFilter(atoms), trajectory='opt.Ag.traj') opt.run(fmax=0.05) print('\n') print(f'Original_Lattice: {original_lattice}') print(f'Final Lattice: {atoms.cell.get_bravais_lattice()}') print(f'Final Cell Volume: {atoms.get_volume():.3f}AA^3') .. rst-class:: sphx-glr-script-out .. code-block:: none Step Time Energy fmax BFGS: 0 02:28:22 0.001584 0.198735 BFGS: 1 02:28:22 0.000255 0.113074 BFGS: 2 02:28:22 -0.000366 0.003126 Original_Lattice: FCC(a=4.09) Final Lattice: FCC(a=4.06315) Final Cell Volume: 16.770AA^3 .. _sphx_glr_download_examples_generated_tutorials_bulk.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: bulk.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: bulk.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: bulk.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_