0
$\begingroup$

I'm looking for a fast way to generate a crystal lattice. I know you can do a bunch of linear combinations of lattice vectors, but this can be costly. Here's my current code for doing so in Python, for a FCC lattice:

import numpy as np

a = 5.472 # lattice parameter

# Lattice vector matrix
A = a*np.array([[0.0, 0.5, 0.5],
                   [0.5, 0.0, 0.5],
                   [0.5, 0.5, 0.0]] )

# Basis in direct coordinates
basisd = np.array ( [[0.0, 0.0, 0.0]] )

# Define the orthohombic box
bx = np.array([ [a, 0., 0.],
                [0., a, 0.],
                [0., 0., a] ] )

# Define the supercell dimensions
supercell = np.array([3., 3., 3.])

# Define maximum scope of n values to check with r= n1*a1 + n2*a2 + n3*a3 (use a larger number for a larger lattice)
nvals = 10;

# Convert to cartesian coordinates
basisc = np.matmul(basisd, A)
bx = np.multiply(bx, supercell)
bx_inv = np.linalg.inv(bx)

a1 = A[0]
a2 = A[1]
a3 = A[2]

basis_shape = np.shape(basisc)
batoms = basis_shape[0]
pos = []
nums = []
for b in range(0,batoms):
    for n1 in range(-nvals, nvals+1):
        for n2 in range(-nvals, nvals+1):
            for n3 in range(-nvals, nvals+1):
                r = n1*a1 + n2*a2 + n3*a3 + basisc[b]

                if 0 <= r[0] and r[0] < bx[0][0] and 0 <= r[1] and r[1] < bx[1][1] and 0 <= r[2] and r[2] < bx[2][2]:

                    pos.append(r)

I had to throw an if-statement in that loop to ensure that I'm only capturing positions that are inside the box, since I'm sampling many positions outside the box. The final array pos has all the lattice positions in the box. This code is fine for a small 3X3X3 unit cell, but it is incredibly slow for larger cells, where I have to sample more lattice vector linear combinations.

Instead of looping over many integers $n_1, n_2, n_3$ to calculate positions $\textbf{r} = n_1\textbf{a}_1+n_2\textbf{a}_2 + n_3\textbf{a}_3$, is there a faster way to do this?

$\endgroup$
4
  • $\begingroup$ This is a question about coding, not physics. Computational Science SE would be a better place to ask such questions. $\endgroup$ May 11, 2018 at 19:50
  • $\begingroup$ You may also want to look at Code Review for some Python-oriented pointers. $\endgroup$
    – Kyle Kanos
    May 11, 2018 at 19:56
  • $\begingroup$ I know people use Python for coding but honestly it's dreadfully inefficient when you have to use Numpy. At least try PyPy or better yet Java or C++. $\endgroup$ May 11, 2018 at 20:39
  • $\begingroup$ @StephenG That nested Python loop is the killer. Generally, if you're using explicit Python loops on Numpy arrays you're doing it wrong, although sometimes such loops are unavoidable. Numpy is fine when it can handle the looping itself at compiled speed, and it's performing arithmetic with machine numeric types. FWIW, converting machine numbers to Python objects is relatively slow, extracting the machine numeric value from an object is quick. $\endgroup$
    – PM 2Ring
    May 12, 2018 at 21:46

1 Answer 1

1
$\begingroup$

You don't seem to be using the symmetry of the crystal to your advantage. In fact you could get away to simply compute coordinates for a single $n_1n_2$ plane, and then make N copies along the $n_3$ direction while adding $\bf{a}_3$ in every iteration to the whole plane

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.