...
 
Commits (2742)
......@@ -7,6 +7,7 @@ MANIFEST
*.traj
/build
/ase.egg-info
/ase/db/static/jsmol
#editor backup files
*.py~
......
services:
- postgres:latest
- mysql:latest
- mariadb:latest
python_2_7_tests:
image: python:2.7
variables:
POSTGRES_DB: testase
POSTGRES_USER: ase
POSTGRES_PASSWORD: "ase"
MYSQL_DATABASE: testase_mysql
MYSQL_ROOT_PASSWORD: ase
# Check oldest supported Python with oldest supported libraries.
# Does not install any optional libraries except matplotlib.
# This test does not use the --strict flag because things like
# deprecation warnings may be rampant, yet these are only important
# for the future, not the past.
#
# Py3.5 wheels exist from scipy-0.16.x, July 2015. We choose scipy-0.16.1.
# Consistent with numpy-1.10.x (October 2015). We choose 1.10.4.
python_3_oldlibs_tests:
image: python:3.5-slim
script:
# manually install ase deps
- pip install numpy scipy matplotlib flask
# extra packages for testing
- pip install pyflakes netCDF4
# using "install from source" instructions
- export PATH=$PATH:$CI_PROJECT_DIR/bin
- echo $PATH
- export PYTHONPATH=$CI_PROJECT_DIR
- echo $PYTHONPATH
# tests
- pip install numpy==1.10.4 scipy==0.16.1 matplotlib==2.0.0
- pip install --no-deps .
- python --version
- ase info
- ase test
# Check newest Python with all the standard dependencies at newest versions.
#
# We currently use scipy 1.2.1 because of failures with 1.3.0 on py3.7.
# TODO: Make things work with 1.3.0
python_3_tests:
image: python:3
image: python:3-slim
script:
# manually install extra deps
#- pip install netCDF4 pyflakes
- pip install pyflakes
# 'pip install' pulls in all required packages
- pip install flake8 psycopg2-binary netCDF4 pymysql cryptography
- pip install flake8 psycopg2-binary netCDF4
- pip install scipy==1.2.1 # TODO: Delete me
- pip install .
# Run tests
- python --version
- ase info
- ase test
# pyflakes code check (Python 3 version)
- ase test --strict
- cd $CI_PROJECT_DIR
- pyflakes ase doc
# Currently (2019-07-26) facing a problem where PATH does not include
# /opt/conda/bin. So we manually update the PATH. This could be a temporary
# issue with the gitlab runners.
#
# Same scipy versioning issue as for the ordinary tests.
conda_tests:
image: continuumio/anaconda3
script:
# manually install extra deps
- apt-get update
- apt-get install -yq libgl1-mesa-glx
- conda install -yq pip wheel numpy scipy pyflakes matplotlib flask
# 'pip install' pulls in all required packages
- echo $PATH
- export PATH=/opt/conda/bin:$PATH
- conda install -yq pip wheel numpy scipy==1.2.1 matplotlib flask
- pip install .
# Run tests
- python --version
- ase info
- ase test
- ase test --strict
docs_test:
image: python:3
script:
# 'pip install' pulls in all required packages (include docs)
- pip install .[docs]
# build the docs
- ase info
- which sphinx-build
- cd $CI_PROJECT_DIR/doc
......
......@@ -70,7 +70,7 @@ Jens Jørgen Mortensen <jensj@fysik.dtu.dk> Jens Jorgen Mortensen <jensj@surt.fy
Jens Jørgen Mortensen <jensj@fysik.dtu.dk> Jens Jørgen Mortensen <jjmo@dtu.dk>
Jesper Friis <jesper.friis@sintef.no> jesperf <jesperf@5af997e4-4c3d-0410-961e-ac0512e437d9>
Jess Wellendorff <mail@example.com> jesswe <jesswe@5af997e4-4c3d-0410-961e-ac0512e437d9>
jg <jag@psibox.dk>
Jakob Gath <jag@psibox.dk> jg <jag@psibox.dk>
Jin Chang <jchang@dtu.dk>
Jingzhe Chen <mail@example.com> jingzhe <jingzhe@5af997e4-4c3d-0410-961e-ac0512e437d9>
John Kitchin <jkitchin@cmu.edu> jk7683@kit.edu <jk7683@uc1n996.localdomain>
......@@ -103,7 +103,6 @@ Maja Lenz <lenz@fhi-berlin.mpg.de> Maja-Olivia Lenz <lenz@fhi-berlin.mpg.de>
Marc Barbry <marc.barbarosa@gmail.com> marc barbry <marc.barbarosa@gmail.com>
Marc Barbry <marc.barbarosa@gmail.com> marc barbry <marc.barbry@mailoo.org>
Marko Melander <marmela@dtu.dk>
markus <markus@5af997e4-4c3d-0410-961e-ac0512e437d9>
Martin Hangaard Hansen <mhah@surt.fysik.dtu.dk> mhah <mhah@5af997e4-4c3d-0410-961e-ac0512e437d9>
Mathias Ljungberg <mathias.ljungberg@gmail.com> Mathias Ljungberg <mathias@captiveportal-10-100-14-106.local.su.se>
Mathias Ljungberg <mathias.ljungberg@gmail.com> Mathias Ljungberg <mathias@captiveportal-10-100-17-25.local.su.se>
......@@ -127,7 +126,7 @@ Morten Gjerding <mogje@fysik.dtu.dk> Morten Gjerding <mogje@fysik.dtu.dk>
Morten Nagel <morten.nagel@helsinki.fi>
Poul Georg Moses <mail@example.com> moses <moses@5af997e4-4c3d-0410-961e-ac0512e437d9>
Marco Vanin <mail@example.com> mvanin <mvanin@5af997e4-4c3d-0410-961e-ac0512e437d9>
mvdb <maxime.cp.vandenbossche@gmail.com>
Maxime van den Bossche <maxime.cp.vandenbossche@gmail.com> mvdb <maxime.cp.vandenbossche@gmail.com>
Nicki Frank Hinsche <nickih@surt.fysik.dtu.dk>
Ole Schütt <ole.schuett@empa.ch> Ole Schuett <ole.schuett@empa.ch>
Ole Schütt <ole.schuett@empa.ch> Ole Schütt <oschuett@cp2k.org>
......@@ -135,6 +134,7 @@ Ole Schütt <ole.schuett@empa.ch> oschuett <oschuett@5af997e4-4c3d-0410-961e-ac0
Oliver Brügner <oliver-bruegner@web.de>
Otto Kohulak <kohulak@Maxwell00> Otto Kohulák <pravod@gmail.com>
Markus Kaukonen <markus.kaukonen@iki.fi> paapu68 <markus.kaukonen@iki.fi>
Markus Kaukonen <markus.kaukonen@iki.fi> markus <markus@5af997e4-4c3d-0410-961e-ac0512e437d9>
Paul C. Jennings <jennings.p.c@gmail.com> Paul C. Jennings <pcje@surt.fysik.dtu.dk>
Paul Erhart <erhart@chalmers.se>
Paweł T. Jochym <pawel.jochym@ifj.edu.pl> Paweł T. Jochym <jochym@wolf.ifj.edu.pl>
......
......@@ -2,7 +2,6 @@ include MANIFEST.in
include COPYING* LICENSE README.rst CONTRIBUTING.rst CHANGELOG.rst
include requirements.txt
include bin/ase
include bin/ase3
include ase/spacegroup/spacegroup.dat
include ase/collections/*.json
include ase/db/static/*
......
......@@ -10,13 +10,16 @@ Webpage: http://wiki.fysik.dtu.dk/ase
Requirements
------------
* Python_ 2.7, 3.4-3.6
* Python_ 3.5 or later
* NumPy_ (base N-dimensional array package)
* SciPy_ (library for scientific computing)
Optional:
* SciPy_ (library for scientific computing)
* For ASE's GUI: Matplotlib_ (2D Plotting)
* tkinter (for ase.gui)
* Flask (for ase.db web-interface)
Installation
......@@ -68,6 +71,8 @@ BFGS: 3 19:10:51 -31.492848 0.0023
>>> h2.get_potential_energy() # ASE's units are eV and Ang
-31.492847800329216
This example requires NWChem to be installed.
::
$ ase gui h2.traj
......
......@@ -3,23 +3,28 @@
"""Atomic Simulation Environment."""
import sys
from distutils.version import LooseVersion
import numpy as np
from ase.atom import Atom
from ase.atoms import Atoms
if sys.version_info[0] == 2:
raise ImportError('ASE requires Python3. This is Python2.')
if LooseVersion(np.__version__) < '1.9':
raise ImportError(
'ASE needs NumPy-1.9.0 or later. You have: %s' % np.__version__)
__all__ = ['Atoms', 'Atom']
__version__ = '3.16.3b1'
__version__ = '3.19.0b1'
from ase.atom import Atom
from ase.atoms import Atoms
# import ase.parallel early to avoid circular import problems when
# ase.parallel does "from gpaw.mpi import world":
import ase.parallel # noqa
ase.parallel # silence pyflakes
if LooseVersion(np.__version__) < '1.9':
# Make isinstance(x, numbers.Integral) work also for np.intxx:
import numbers
numbers.Integral.register(np.integer)
del numbers
This diff is collapsed.
from ase.cli.build import main
main()
......@@ -2,7 +2,8 @@ from optparse import OptionParser
import numpy as np
from ase.atoms import Atoms, string2symbols
from ase.atoms import Atoms
from ase.symbols import string2symbols
from ase.build import (molecule, fcc111, hcp0001, bcc110, bcc100, diamond111,
add_adsorbate)
from ase.data import reference_states, atomic_numbers, covalent_radii
......
from __future__ import division
from math import sqrt
from ase.atoms import Atoms, string2symbols
from ase.atoms import Atoms
from ase.symbols import string2symbols
from ase.data import reference_states, atomic_numbers, chemical_symbols
from ase.utils import plural
......@@ -46,7 +47,7 @@ def bulk(name, crystalstructure=None, a=None, c=None, covera=None, u=None,
xref = ref['symmetry']
structures = {'sc': 1, 'fcc': 1, 'bcc': 1, 'hcp': 1, 'diamond': 1,
'zincblende': 2, 'rocksalt':2, 'cesiumchloride':2,
'zincblende': 2, 'rocksalt': 2, 'cesiumchloride': 2,
'fluorite': 3, 'wurtzite': 2}
if crystalstructure is None:
......
......@@ -19,6 +19,7 @@ def surface(lattice, indices, layers, vacuum=None, tol=1e-10):
Number of equivalent layers of the slab.
vacuum: float
Amount of vacuum added on both sides of the slab.
"""
indices = np.asarray(indices)
......@@ -31,6 +32,7 @@ def surface(lattice, indices, layers, vacuum=None, tol=1e-10):
h, k, l = indices
h0, k0, l0 = (indices == 0)
if h0 and k0 or h0 and l0 or k0 and l0: # if two indices are zero
if not h0:
c1, c2, c3 = [(0, 1, 0), (0, 0, 1), (1, 0, 0)]
......
......@@ -17,7 +17,8 @@ def graphene_nanoribbon(n, m, type='zigzag', saturated=False, C_H=1.09,
Parameters:
n: int
The width of the nanoribbon.
The width of the nanoribbon. For armchair nanoribbons, this
n may be half-integer to repeat by half a cell.
m: int
The length of the nanoribbon.
type: str
......@@ -39,6 +40,11 @@ def graphene_nanoribbon(n, m, type='zigzag', saturated=False, C_H=1.09,
If true, make an infinite sheet instead of a ribbon (default: False)
"""
if m % 1 != 0:
raise ValueError('m must be integer')
if type == 'zigzag' and n % 1 != 0:
raise ValueError('n must be an integer for zigzag ribbons')
b = sqrt(3) * C_C / 4
arm_unit = Atoms(main_element + '4',
pbc=(1, 0, 1),
......@@ -47,6 +53,11 @@ def graphene_nanoribbon(n, m, type='zigzag', saturated=False, C_H=1.09,
[b * 2, 0, C_C / 2.],
[b * 2, 0, 3 * C_C / 2.],
[0, 0, 2 * C_C]]
arm_unit_half = Atoms(main_element + '2',
pbc=(1, 0, 1),
cell=[2 * b, 0, 3 * C_C])
arm_unit_half.positions = [[b * 2, 0, C_C / 2.],
[b * 2, 0, 3 * C_C / 2.]]
zz_unit = Atoms(main_element + '2',
pbc=(1, 0, 1),
cell=[3 * C_C / 2.0, 0, b * 4])
......@@ -88,30 +99,50 @@ def graphene_nanoribbon(n, m, type='zigzag', saturated=False, C_H=1.09,
atoms.cell = [n * 3 * C_C / 2, 0, m * 4 * b]
elif type == 'armchair':
for i in range(n):
n *= 2
n_int = int(round(n))
if abs(n_int - n) > 1e-10:
raise ValueError(
'The argument n has to be half-integer for armchair ribbons.')
n = n_int
for i in range(n // 2):
layer = arm_unit.repeat((1, 1, m))
layer.positions[:, 0] -= 4 * b * i
atoms += layer
if n % 2:
layer = arm_unit_half.repeat((1, 1, m))
layer.positions[:, 0] -= 4 * b * (n // 2)
atoms += layer
xmin = atoms.positions[-1, 0]
if saturated:
arm_right_saturation = Atoms(saturate_element + '2', pbc=(1, 0, 1),
cell=[4 * b, 0, 3 * C_C])
arm_right_saturation.positions = [
[- sqrt(3) / 2 * C_H, 0, C_H * 0.5],
[- sqrt(3) / 2 * C_H, 0, 2 * C_C - C_H * 0.5]]
if n % 2:
arm_right_saturation = Atoms(saturate_element + '2',
pbc=(1, 0, 1),
cell=[2 * b, 0, 3 * C_C])
arm_right_saturation.positions = [
[- sqrt(3) / 2 * C_H, 0, C_C / 2 - C_H * 0.5],
[- sqrt(3) / 2 * C_H, 0, 3 * C_C / 2.0 + C_H * 0.5]]
else:
arm_right_saturation = Atoms(saturate_element + '2',
pbc=(1, 0, 1),
cell=[4 * b, 0, 3 * C_C])
arm_right_saturation.positions = [
[- sqrt(3) / 2 * C_H, 0, C_H * 0.5],
[- sqrt(3) / 2 * C_H, 0, 2 * C_C - C_H * 0.5]]
arm_left_saturation = Atoms(saturate_element + '2', pbc=(1, 0, 1),
cell=[4 * b, 0, 3 * C_C])
arm_left_saturation.positions = [
[b * 2 + sqrt(3) / 2 * C_H, 0, C_C / 2 - C_H * 0.5],
[b * 2 + sqrt(3) / 2 * C_H, 0, 3 * C_C / 2.0 + C_H * 0.5]]
arm_right_saturation.positions[:, 0] -= 4 * b * (n - 1)
arm_right_saturation.positions[:, 0] -= 4 * b * (n / 2.0 - 1)
atoms += arm_right_saturation.repeat((1, 1, m))
atoms += arm_left_saturation.repeat((1, 1, m))
atoms.cell = [b * 4 * n, 0, 3 * C_C * m]
atoms.cell = [b * 4 * n / 2.0, 0, 3 * C_C * m]
atoms.set_pbc([sheet, False, True])
......@@ -120,7 +151,7 @@ def graphene_nanoribbon(n, m, type='zigzag', saturated=False, C_H=1.09,
atoms.positions[:, 0] -= xmin
if not sheet:
atoms.cell[0] = 0.0
if vacuum:
if vacuum is not None:
atoms.center(vacuum, axis=1)
if not sheet:
atoms.center(vacuum, axis=0)
......
......@@ -2,8 +2,14 @@
import numpy as np
from ase import Atoms
def get_deviation_from_optimal_cell_shape(cell, target_shape='sc', norm=None):
class SupercellError(Exception):
"""Use if construction of supercell fails"""
def get_deviation_from_optimal_cell_shape(cell, target_shape="sc", norm=None):
"""
Calculates the deviation of the given cell metric from the ideal
cell metric defining a certain shape. Specifically, the function
......@@ -28,21 +34,25 @@ def get_deviation_from_optimal_cell_shape(cell, target_shape='sc', norm=None):
"""
if target_shape in ['sc', 'simple-cubic']:
if target_shape in ["sc", "simple-cubic"]:
target_metric = np.eye(3)
elif target_shape in ['fcc', 'face-centered cubic']:
target_metric = 0.5 * np.array([[0, 1, 1],
[1, 0, 1],
[1, 1, 0]])
elif target_shape in ["fcc", "face-centered cubic"]:
target_metric = 0.5 * np.array([[0, 1, 1], [1, 0, 1], [1, 1, 0]])
if not norm:
norm = (np.linalg.det(cell) /
np.linalg.det(target_metric))**(-1.0 / 3)
norm = (np.linalg.det(cell) / np.linalg.det(target_metric)) ** (
-1.0 / 3
)
return np.linalg.norm(norm * cell - target_metric)
def find_optimal_cell_shape(cell, target_size, target_shape,
lower_limit=-2, upper_limit=2,
verbose=False):
def find_optimal_cell_shape(
cell,
target_size,
target_shape,
lower_limit=-2,
upper_limit=2,
verbose=False,
):
"""Returns the transformation matrix that produces a supercell
corresponding to *target_size* unit cells with metric *cell* that
most closely approximates the shape defined by *target_shape*.
......@@ -67,35 +77,37 @@ def find_optimal_cell_shape(cell, target_size, target_shape,
"""
# Set up target metric
if target_shape in ['sc', 'simple-cubic']:
if target_shape in ["sc", "simple-cubic"]:
target_metric = np.eye(3)
elif target_shape in ['fcc', 'face-centered cubic']:
target_metric = 0.5 * np.array([[0, 1, 1],
[1, 0, 1],
[1, 1, 0]], dtype=float)
elif target_shape in ["fcc", "face-centered cubic"]:
target_metric = 0.5 * np.array(
[[0, 1, 1], [1, 0, 1], [1, 1, 0]], dtype=float
)
if verbose:
print('target metric (h_target):')
print("target metric (h_target):")
print(target_metric)
# Normalize cell metric to reduce computation time during looping
norm = (target_size * np.linalg.det(cell) /
np.linalg.det(target_metric))**(-1.0 / 3)
norm = (
target_size * np.linalg.det(cell) / np.linalg.det(target_metric)
) ** (-1.0 / 3)
norm_cell = norm * cell
if verbose:
print('normalization factor (Q): %g' % norm)
print("normalization factor (Q): %g" % norm)
# Approximate initial P matrix
ideal_P = np.dot(target_metric, np.linalg.inv(norm_cell))
if verbose:
print('idealized transformation matrix:')
print("idealized transformation matrix:")
print(ideal_P)
starting_P = np.array(np.around(ideal_P, 0), dtype=int)
if verbose:
print('closest integer transformation matrix (P_0):')
print("closest integer transformation matrix (P_0):")
print(starting_P)
# Prepare run.
from itertools import product
best_score = 1e6
optimal_P = None
for dP in product(range(lower_limit, upper_limit + 1), repeat=9):
......@@ -104,28 +116,31 @@ def find_optimal_cell_shape(cell, target_size, target_shape,
if int(np.around(np.linalg.det(P), 0)) != target_size:
continue
score = get_deviation_from_optimal_cell_shape(
np.dot(P, norm_cell), target_shape=target_shape, norm=1.0)
np.dot(P, norm_cell), target_shape=target_shape, norm=1.0
)
if score < best_score:
best_score = score
optimal_P = P
if optimal_P is None:
print('Failed to find a transformation matrix.')
print("Failed to find a transformation matrix.")
return None
# Finalize.
if verbose:
print('smallest score (|Q P h_p - h_target|_2): %f' % best_score)
print('optimal transformation matrix (P_opt):')
print("smallest score (|Q P h_p - h_target|_2): %f" % best_score)
print("optimal transformation matrix (P_opt):")
print(optimal_P)
print('supercell metric:')
print("supercell metric:")
print(np.round(np.dot(optimal_P, cell), 4))
print('determinant of optimal transformation matrix: %g' %
np.linalg.det(optimal_P))
print(
"determinant of optimal transformation matrix: %g"
% np.linalg.det(optimal_P)
)
return optimal_P
def make_supercell(prim, P):
def make_supercell(prim, P, wrap=True, tol=1e-5):
"""Generate a supercell by applying a general transformation (*P*) to
the input configuration (*prim*).
......@@ -135,16 +150,92 @@ def make_supercell(prim, P):
configuraton `\mathbf{h}_p` by `\mathbf{P h}_p =
\mathbf{h}`.
Internally this function uses the :func:`~ase.build.cut` function.
Parameters:
prim: ASE Atoms object
Input configuration.
P: 3x3 integer matrix
Transformation matrix `\mathbf{P}`.
wrap: bool
wrap in the end
tol: float
tolerance for wrapping
"""
supercell_matrix = P
supercell = clean_matrix(supercell_matrix @ prim.cell)
# cartesian lattice points
lattice_points_frac = lattice_points_in_supercell(supercell_matrix)
lattice_points = np.dot(lattice_points_frac, supercell)
superatoms = Atoms(cell=supercell, pbc=prim.pbc)
for lp in lattice_points:
shifted_atoms = prim.copy()
shifted_atoms.positions += lp
superatoms.extend(shifted_atoms)
# check number of atoms is correct
n_target = int(np.round(np.linalg.det(supercell_matrix) * len(prim)))
if n_target != len(superatoms):
msg = "Number of atoms in supercell: {}, expected: {}".format(
n_target, len(superatoms)
)
raise SupercellError(msg)
if wrap:
superatoms.wrap(eps=tol)
return superatoms
def lattice_points_in_supercell(supercell_matrix):
"""Find all lattice points contained in a supercell.
Adapted from pymatgen, which is available under MIT license:
The MIT License (MIT) Copyright (c) 2011-2012 MIT & The Regents of the
University of California, through Lawrence Berkeley National Laboratory
"""
from ase.build import cut
return cut(prim, P[0], P[1], P[2])
diagonals = np.array(
[
[0, 0, 0],
[0, 0, 1],
[0, 1, 0],
[0, 1, 1],
[1, 0, 0],
[1, 0, 1],
[1, 1, 0],
[1, 1, 1],
]
)
d_points = np.dot(diagonals, supercell_matrix)
mins = np.min(d_points, axis=0)
maxes = np.max(d_points, axis=0) + 1
ar = np.arange(mins[0], maxes[0])[:, None] * np.array([1, 0, 0])[None, :]
br = np.arange(mins[1], maxes[1])[:, None] * np.array([0, 1, 0])[None, :]
cr = np.arange(mins[2], maxes[2])[:, None] * np.array([0, 0, 1])[None, :]
all_points = ar[:, None, None] + br[None, :, None] + cr[None, None, :]
all_points = all_points.reshape((-1, 3))
frac_points = np.dot(all_points, np.linalg.inv(supercell_matrix))
tvects = frac_points[
np.all(frac_points < 1 - 1e-10, axis=1)
& np.all(frac_points >= -1e-10, axis=1)
]
assert len(tvects) == round(abs(np.linalg.det(supercell_matrix)))
return tvects
def clean_matrix(matrix, eps=1e-12):
""" clean from small values"""
matrix = np.array(matrix)
for ij in np.ndindex(matrix.shape):
if abs(matrix[ij]) < eps:
matrix[ij] = 0
return matrix
......@@ -18,17 +18,21 @@ from ase.lattice.cubic import FaceCenteredCubic
from ase.utils import basestring
def fcc100(symbol, size, a=None, vacuum=None, orthogonal=True):
def fcc100(symbol, size, a=None, vacuum=None, orthogonal=True,
periodic=False):
"""FCC(100) surface.
Supported special adsorption sites: 'ontop', 'bridge', 'hollow'."""
if not orthogonal:
raise NotImplementedError("Can't do non-orthogonal cell yet!")
return _surface(symbol, 'fcc', '100', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'fcc', '100', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def fcc110(symbol, size, a=None, vacuum=None, orthogonal=True):
def fcc110(symbol, size, a=None, vacuum=None, orthogonal=True,
periodic=False):
"""FCC(110) surface.
Supported special adsorption sites: 'ontop', 'longbridge',
......@@ -36,20 +40,27 @@ def fcc110(symbol, size, a=None, vacuum=None, orthogonal=True):
if not orthogonal:
raise NotImplementedError("Can't do non-orthogonal cell yet!")
return _surface(symbol, 'fcc', '110', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'fcc', '110', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def bcc100(symbol, size, a=None, vacuum=None, orthogonal=True):
def bcc100(symbol, size, a=None, vacuum=None, orthogonal=True,
periodic=False):
"""BCC(100) surface.
Supported special adsorption sites: 'ontop', 'bridge', 'hollow'."""
if not orthogonal:
raise NotImplementedError("Can't do non-orthogonal cell yet!")
return _surface(symbol, 'bcc', '100', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'bcc', '100', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def bcc110(symbol, size, a=None, vacuum=None, orthogonal=False):
def bcc110(symbol, size, a=None, vacuum=None, orthogonal=False,
periodic=False):
"""BCC(110) surface.
Supported special adsorption sites: 'ontop', 'longbridge',
......@@ -57,40 +68,52 @@ def bcc110(symbol, size, a=None, vacuum=None, orthogonal=False):
Use *orthogonal=True* to get an orthogonal unit cell - works only
for size=(i,j,k) with j even."""
return _surface(symbol, 'bcc', '110', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'bcc', '110', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def bcc111(symbol, size, a=None, vacuum=None, orthogonal=False):
def bcc111(symbol, size, a=None, vacuum=None, orthogonal=False,
periodic=False):
"""BCC(111) surface.
Supported special adsorption sites: 'ontop'.
Use *orthogonal=True* to get an orthogonal unit cell - works only
for size=(i,j,k) with j even."""
return _surface(symbol, 'bcc', '111', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'bcc', '111', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def fcc111(symbol, size, a=None, vacuum=None, orthogonal=False):
def fcc111(symbol, size, a=None, vacuum=None, orthogonal=False,
periodic=False):
"""FCC(111) surface.
Supported special adsorption sites: 'ontop', 'bridge', 'fcc' and 'hcp'.
Use *orthogonal=True* to get an orthogonal unit cell - works only
for size=(i,j,k) with j even."""
return _surface(symbol, 'fcc', '111', size, a, None, vacuum, orthogonal)
return _surface(symbol, 'fcc', '111', size, a, None, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def hcp0001(symbol, size, a=None, c=None, vacuum=None, orthogonal=False):
def hcp0001(symbol, size, a=None, c=None, vacuum=None, orthogonal=False,
periodic=False):
"""HCP(0001) surface.
Supported special adsorption sites: 'ontop', 'bridge', 'fcc' and 'hcp'.
Use *orthogonal=True* to get an orthogonal unit cell - works only
for size=(i,j,k) with j even."""
return _surface(symbol, 'hcp', '0001', size, a, c, vacuum, orthogonal)
return _surface(symbol, 'hcp', '0001', size, a, c, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def hcp10m10(symbol, size, a=None, c=None, vacuum=None, orthogonal=True):
def hcp10m10(symbol, size, a=None, c=None, vacuum=None, orthogonal=True,
periodic=False):
"""HCP(10m10) surface.
Supported special adsorption sites: 'ontop'.
......@@ -99,10 +122,13 @@ def hcp10m10(symbol, size, a=None, c=None, vacuum=None, orthogonal=True):
if not orthogonal:
raise NotImplementedError("Can't do non-orthogonal cell yet!")
return _surface(symbol, 'hcp', '10m10', size, a, c, vacuum, orthogonal)
return _surface(symbol, 'hcp', '10m10', size, a, c, vacuum,
periodic=periodic,
orthogonal=orthogonal)
def diamond100(symbol, size, a=None, vacuum=None, orthogonal=True):
def diamond100(symbol, size, a=None, vacuum=None, orthogonal=True,
periodic=False):
"""DIAMOND(100) surface.
Supported special adsorption sites: 'ontop'."""
......@@ -110,10 +136,12 @@ def diamond100(symbol, size, a=None, vacuum=None, orthogonal=True):
raise NotImplementedError("Can't do non-orthogonal cell yet!")
return _surface(symbol, 'diamond', '100', size, a, None, vacuum,
orthogonal)
periodic=periodic,
orthogonal=orthogonal)
def diamond111(symbol, size, a=None, vacuum=None, orthogonal=False):
def diamond111(symbol, size, a=None, vacuum=None, orthogonal=False,
periodic=False):
"""DIAMOND(111) surface.
Supported special adsorption sites: 'ontop'."""
......@@ -121,7 +149,8 @@ def diamond111(symbol, size, a=None, vacuum=None, orthogonal=False):
if orthogonal:
raise NotImplementedError("Can't do orthogonal cell yet!")
return _surface(symbol, 'diamond', '111', size, a, None, vacuum,
orthogonal)
periodic=periodic,
orthogonal=orthogonal)
def add_adsorbate(slab, adsorbate, height, position=(0, 0), offset=None,
......@@ -242,7 +271,8 @@ def add_vacuum(atoms, vacuum):
atoms.set_cell(uc)
def _surface(symbol, structure, face, size, a, c, vacuum, orthogonal=True):
def _surface(symbol, structure, face, size, a, c, vacuum, periodic,
orthogonal=True):
"""Function to build often used surfaces.
Don't call this function directly - use fcc100, fcc110, bcc111, ..."""
......@@ -274,7 +304,7 @@ def _surface(symbol, structure, face, size, a, c, vacuum, orthogonal=True):
slab = Atoms(numbers,
tags=tags.ravel(),
pbc=(True, True, False),
pbc=(True, True, periodic),
cell=size)
surface_cell = None
......@@ -379,9 +409,10 @@ def _surface(symbol, structure, face, size, a, c, vacuum, orthogonal=True):
cell = np.diag(cell)
slab.set_positions(positions.reshape((-1, 3)))
slab.set_cell([a * v * n for v, n in zip(cell, size)], scale_atoms=True)
slab.cell[2] = 0.0
if not periodic:
slab.cell[2] = 0.0
if vacuum is not None:
slab.center(vacuum, axis=2)
......@@ -391,7 +422,6 @@ def _surface(symbol, structure, face, size, a, c, vacuum, orthogonal=True):
slab.info['adsorbate_info']['cell'] = surface_cell
slab.info['adsorbate_info']['sites'] = sites
return slab
......@@ -475,3 +505,12 @@ def mx2(formula='MoS2', kind='2H', a=3.18, thickness=3.19,
atoms = atoms.repeat(size)
return atoms
def _all_surface_functions():
# Convenient for debugging.
d = {}
for func in [fcc100, fcc110, bcc100, bcc110, bcc111, fcc111, hcp0001,
hcp10m10, diamond100, diamond111, fcc111, mx2]:
d[func.__name__] = func
return d
import numpy as np
from ase.build.general_surface import surface
from ase.geometry import get_layers
from ase.symbols import string2symbols
def surfaces_with_termination(lattice, indices, layers, vacuum=None, tol=1e-10,
termination=None, return_all=False, verbose=False):
"""Create surface from a given lattice and Miller indices with a given
termination
Parameters
==========
lattice: Atoms object or str
Bulk lattice structure of alloy or pure metal. Note that the
unit-cell must be the conventional cell - not the primitive cell.
One can also give the chemical symbol as a string, in which case the
correct bulk lattice will be generated automatically.
indices: sequence of three int
Surface normal in Miller indices (h,k,l).
layers: int
Number of equivalent layers of the slab. (not the same as the layers
you choose from for terminations)
vacuum: float
Amount of vacuum added on both sides of the slab.
termination: str
the atoms you wish to be in the top layer. There may be many such
terminations, this function returns all terminations with the same
atomic composition.
e.g. 'O' will return oxygen terminated surfaces.
e.g.'TiO' will return surfaces terminated with layers containing both O
and Ti
Returns:
return_surfs: List
a list of surfaces that match the specifications given
"""
lats = translate_lattice(lattice, indices)
return_surfs = []
check = []
check2 = []
for item in lats:
too_similar = False
surf = surface(item, indices, layers, vacuum=vacuum, tol=tol)
surf.wrap(pbc = [True] * 3) # standardize slabs
positions = surf.get_scaled_positions().flatten()
for i, value in enumerate(positions):
if value >= 1 - tol: # move things closer to zero within tol
positions[i] -= 1
surf.set_scaled_positions(np.reshape(positions, (len(surf), 3)))
#rep = find_z_layers(surf)
z_layers, hs = get_layers(surf, (0, 0, 1)) # just z layers matter
# get the indicies of the atoms in the highest layer
top_layer = [i for i, val in enumerate(z_layers == max(z_layers)) if val]
if termination is not None:
comp = [surf.get_chemical_symbols()[a] for a in top_layer]
term = string2symbols(termination)
# list atoms in top layer and not in requested termination
check = [a for a in comp if a not in term]
# list of atoms in requested termination and not in top layer
check2 = [a for a in term if a not in comp]
if len(return_surfs) > 0:
pos_diff = [a.get_positions() - surf.get_positions()
for a in return_surfs]
for i, su in enumerate(pos_diff):
similarity_test = su.flatten() < tol * 1000
if similarity_test.all():
# checks if surface is too similar to another surface
too_similar = True
if too_similar:
continue
if return_all is True:
pass
elif check != [] or check2 != []:
continue
return_surfs.append(surf)
return return_surfs
def translate_lattice(lattice, indices, tol=10**-3):
"""translates a bulk unit cell along a normal vector given by the a set of
miller indices to the next symetric position. This is used to control the
termination of the surface in the smart_surface command
Parameters:
==========
lattice: Atoms object
atoms object of the bulk unit cell
indices: 1x3 list,tuple, or numpy array
the miller indices you wish to cut along.
returns:
lattice_list: list of Atoms objects
a list of all the different translations of the unit cell that will
yield different terminations of a surface cut along the miller
indices provided.
"""
lattice_list = []
cell = lattice.get_cell()
pt = [0, 0, 0]
h, k, l = indices