The BCC lattice

The reason why I especially posted this is to show how mathematics can be seamlessly converted into a program. This program generates a periodic lattice. Although the output of this program might seem similar to that shown in a couple of posts prior, this is more closely related to the physics. Methods of translating the mathematical formulation into Python, is shown below.

The Lattice:
A lattice is defined by a set of points given by p(n₁,n₂,n₃) = n₁*a₁ + n₂*a₂ + n₃*a₃.
Here, the three vectors are three unit vectors (not necessarily orthogonal), and the numbers are integer indices to the vectors. The function 'lattice_points(N)' returns a bunch of these integers (n₁,n₂,n₃), such that we may be able to map atoms to these points. 

The Basis:
To each lattice point, a number of atoms are attached, displaced by a certain amount from the lattice point. For example, the BCC lattice has a basis of 2. In the program, the function 'unit_vectors_BCC(distance)' provides the position of the basis atoms relative to a particular point. Notice that the way in which a basis point is defined, one no longer needs to separately define unit vectors. These are defined here. However, they may easily be defined separately for greater modularity. And will probably be necessary in structures that are not cubic. 

The Crystal as a Whole:
For finding the location of each of the position of the atoms within a lattice, one just needs to add all of the basis atoms in each of the lattice sites. This is done in the function 'atom_locations(basis, lattice)'.

Once the atom locations are known, it is just a matter of displaying the spheres.

The result is shown below, for a BCC lattice:

The program is shown below:

######################################

# This program is used for calculating

# the bravis lattice points of a unit

# cell, and displaying the results in

# 3D. 

#

# Note that, since the mathematics of

# the bravis lattices are so uniform,

# it is relatively easy to calculate

# the positions of the positions of 

# all atoms, given a particular basis

######################################

##### Functions for doing the physics ##########

#-----------------------------------------------

# This function iterates over all lattice points

# and produces a list of 3D lattice points

#

# Input: 

#    N - the number of points in each dimension

# Output:

#    a [3 by (N+1)^3] matrix containing all

# the lattice points. Thus, for N = 1, returned

# points are,

# [[0,0,0], [0,0,1], [0,1,0], [0,1,1], [1,0,0], 

#  [1,0,1], [1,1,0], [1,1,1]]

#-----------------------------------------------

def lattice_points(N):

lattice = []

for i in range(N+1):

for j in range(N+1):

for k in range(N+1):

lattice.append([i, j, k])

lattice.append([i,j,k])

return lattice

#----------------------------------------------

# This function returns the unit vectors in the 

# three different directions, for each of the 

# basis atoms. This function is generally so 

# simple that in most situatons, it may be 

# defined by hand. Maybe a class is in order?

#

# Shown below is the unit vectors for a BCC 

# lattice, with lattice size 'distance'.

# Notice that BCC has 2 atoms in its unit cell.

# Thus you only provide only 2 atoms. 

#----------------------------------------------

def unit_vectors_BCC(distance):

a = distance; 

point1 = [a, a, a]

point2 = [a/2, a/2, a/2]

return [point1, point2]

#----------------------------------------------

# This function returns a list of coordinates

# for each basis, given the basis, and the 

# lattice points. 

#----------------------------------------------

def atom_locations(basis, lattice):

# Atom definitions

atom1 = []

atom2 = []

# Basis vectors

a1, a2, a3 = basis[0]

b1, b2, b3 = basis[1]

print a1, a2, a3

print b1, b2, b3

for n1, n2, n3 in lattice:

atom1.append( [n1*a1, n2*a2, n3*a3] )

atom2.append( [n1*a1 + b1, n2*a2 + b2, n3*a3 + b3] )

return [atom1, atom2]

#############################################################

#   At this point, we have all the atom locations. All that 

# remains to do is to draw them

#############################################################

from vtk import *

# create a rendering window and renderer

ren = vtkRenderer()

renWin = vtkRenderWindow()

renWin.AddRenderer(ren)

renWin.SetSize(300,300)

renWin.StereoCapableWindowOn()

iren = vtkRenderWindowInteractor()

iren.SetRenderWindow(renWin)

lattice_size = 10e-10; #meters

atom1_radius = 4e-10;  #meters

atom2_radius = 4e-10;  #meters

latticePoints = lattice_points(2)

a1, a2 = atom_locations(unit_vectors_BCC(lattice_size), latticePoints)

# Create the spheres ...

for a, b in zip(a1, a2):

##### This is for atom 1 ################

# create the object to be drawn

sphere = vtkSphereSource()

sphere.SetRadius(atom1_radius)

sphere.SetThetaResolution(30)

sphere.SetPhiResolution(30)

sphere.SetCenter(a[0], a[1], a[2])

# Convert the sphere into polygons

sphereMapper = vtkPolyDataMapper()

sphereMapper.SetInput(sphere.GetOutput())

#Create an actor for the sphere 

sphereActor = vtkActor()

sphereActor.SetMapper(sphereMapper)

(sphereActor.GetProperty()).SetColor(1,0,0)

# assign our actor to the renderer

ren.AddActor(sphereActor)

##### This is for atom 2 ################

# create the object to be drawn

sphere = vtkSphereSource()

sphere.SetRadius(atom2_radius)

sphere.SetThetaResolution(30)

sphere.SetPhiResolution(30)

sphere.SetCenter(b[0], b[1], b[2])

# Convert the sphere into polygons

sphereMapper = vtkPolyDataMapper()

sphereMapper.SetInput(sphere.GetOutput())

#Create an actor for the sphere 

sphereActor = vtkActor()

sphereActor.SetMapper(sphereMapper)

(sphereActor.GetProperty()).SetColor(0,0,1)

(sphereActor.GetProperty()).SetOpacity(1)

# assign our actor to the renderer

ren.AddActor(sphereActor)

# enable user interface interactor

iren.Initialize()

iren.Start()