Dispersion curves for nearest neighbour FCC lattice in the 100 CUBIC direction

IMPORTANT NOTE: The FCC basis was assumed to be a 4-atom cubic basis & the system considered was a cube with periodic boundaries. The dynamical matrices were only calculated for modes in the 100 direction.

Fig 1: The dispersion curves for all 12 modes.

Fig 2: Dispersion for acoustic mode.

Density of states for 2D slices of a 3D harmonic crystal

Fig 1(top) & 2(bottom): These demonstrate that the density of states scales like the 3rd power of frequency.


Fig 3: The 111 patch in this case crosses the system's periodic boundaries (it doesn't self intersect). The DOS still shows the same scaling as Figs 1 & 2.

Computing Green's function for a 2D slice of a 3D crystal

Crystal: FCC lattice with harmonic nearest neighbour interactions. The domain is assumed to be cubic with periodic boundaries.

FCC Basis: 4 atom Cubic FCC basis. The lattice vectors are simply aligned along the cubic axes.

Calculation Scheme: First the 12x12 Hessian matrices are written down as a function of CUBIC lattice separation. Note: since we have only nearest neighbour interactions, we have at most 27 lattice translations whose corresponding hessian matrices have non-zero entries.

The Dynamical matrices, D(k), are then computing for all CUBIC wavevectors that can fit in the simulation box via taking the DFT of Hessian matrices, that is: D(k) = SUM{{R} H(R)*exp(i k.R)}}, where (R} is the set of all possible translations in the lattice (NOTE: R and -R must be separately counted). This DFT was not done through the FFT since the number of non-zero Hessians is much smaller than the number of lattice translations.

Each dynamical matrix was diagonalized to get the 12 polarization vectors & their corresponding eigenvalues {(P_i, L_i) : 1<= i <=12} and the Green's function in fourier space is computed using: G(k) = SUM{(i=1...12) P_i.outerProduct.P_i/L_i }

Finally the real space Green's function can be obtained for the whole crystal through an inverse Fourier Transform: G(R) = SUM{{k} G(k)*exp(i k.R)} where {k} is the set of all fourier modes that can fit in the lattice.

To obtain the Green's function for the slice, I rotated the tensor obtained and discarded all components that had a 'z' piece. Finally I was left with a 2x2 matrix which gave the 2D correlations in a slice (for a particular separation vector) from a 3D crystal.


Difficulties: Due to the fact, that G(k) is the inverse of D(k), we have problems when computing G(k) for long wavelength fourier modes (since their translational and longitudinal polarizations; eigenvalues are very small, hence these D(k) have high condition number). Specifically we can have G(k) & G(-k) differ significantly enough for G(R) to be not purely real. I solved this problem by explicitly setting: G(k)=G(-k).

Density of States - Boson peak

Fig 1: Convergence of DOS of a 111 plane slice of a 3D FCC crystal (20% disorder width). Note the boson peak at around omega=6.9

Convergence of MSD of soft-sphere particles

In the figure above the x-axis is sampling time while the y-axis is the variance of the MSD of 111 plane computed at a particular timestep. The distribution of MSD of the ordered system must go to a delta function as we sample longer. The variance must go like ~ 1/N (where N=number of independent samples). The measured exponent of the power-law decay in the variance of the ordered system's MSD is ~= 1, in accord with expectation. The disordered system's MSD distribution must have an intrinsic variance which is exhibited by the plateauing of the variance.

FCC lattice: planewave projections of full lattice and of 111,100 plane slices

FCC lattice: Cubic Indexing and mapping to Bravais Index, handling PBCs

For constructing neighbour lists, while accounting for PBCs, we have applied a scheme of using 'cubic indexes' for each atom in our FCC lattice.

NOTE: Our simulation cell is cubic, but it could be cuboidal too & this indexing would still be good.

Cubic Index: Each atom is assigned a unique 4-tuple (m,n,p,s) of integers.

The 3-tuple (m,n,p) indexes the cubic Bravais cell. Each cubic cell has 4 atoms which are distinguished by their 's' index according to the following rule: (atom_position) --- ('s' index)
v0: (0,0,0) -- s=0
v1: (1/2,1/2,0) -- s=1
v2: (0,1/2,1/2) -- s=2
v3: (1/2,0,1/2) -- s=3
(v1,v2,v3 are the FCC basis vectors)

Each time a lattice is constructed with cubic indexing, a pickled dict: 'cubicIndexMap.pickle' is specifying a map from (m,n,p,s)-->ID is also written.

A map from the integer Bravais indices (a,b,c), which are integers specified along the primitive FCC basis vectors, can be specified to the cubic indices using only integer math.
F: (a,b,c)-->(m,n,p,s)
(To construct F we write the position of the atom first in Bravais space and then expand the basis vectors into 3D cartesian space, as specified above, noting that EITHER 2 vector components OR none can at most be half-integers - depending upon the value of s).

The map F: (a,b,c)-->(m,n,p,s) is,
m,v[1] = divmod(a+c,2)
n,v[2] = divmod(b+a,2)
p,v[3] = divmod(b+c,2)

s is assigned by checking if 'v' is v0,v1,v2 or v3.
NOTE: divmod () is a built-in python function which returns the quotient & the remainder of a division.

Applying PBCs:
When finding Bravais neighbours, one can quickly map the cubic index and then to the ID (using the pickled map). In case periodic boundaries need to be handled, we merely need to 'mod' cubic indices (m,n,p) with the number of cubic cells along a cartesian axis. NOTE: 's' should be unchanged.
See: the function getCubicIndex() in
/home/ahasan/usr/local/simulations/fcc-3d/disordered/lattice/pickFCCplane.py