3.016 Problem Set #2.Group, 2010

Import data

Import raw data set A, B and C from the web, assign the data sets as aData, bData and cData:

aData=Import["http://pruffle.mit.edu/3.016/Data/atomPositions_A.csv"];

bData=Import["http://pruffle.mit.edu/3.016/Data/atomPositions_B.csv"];

cData=Import["http://pruffle.mit.edu/3.016/Data/atomPositions_C.csv"];

Now checking the dimensions of the imported data:

Dimensions[aData]

Dimensions[bData]

Dimensions[cData]

There are 4905, 4675 and 2581 atoms' 3D position data in the three data sets.

Exploring data sets A, B and C

Let's first visualize the data sets

aplot=ListPointPlot3D[aData,BoxRatios-> {1,1,1}]

bplot=ListPointPlot3D[bData,BoxRatios-> {1,1,1}]

cplot=ListPointPlot3D[cData,BoxRatios-> {1,1,1}]

Trying to rotate each of the above plotted images, one can see that A looks quite random while for both B and C we can see at certain angles something like ordered lattice planes.

Let's visualize as small spheres instead of points... play with the sphere sizes and rotate the images to see if we can find anything interesting

Graphics3D[
GraphicsComplex[
aData,
Table[Sphere[i, 0.15],{i,1,Length[aData]}]
]
]

Graphics3D[
GraphicsComplex[
bData,
Table[Sphere[i, 0.1],{i,1,Length[bData]}]
]
]

Graphics3D[
GraphicsComplex[
cData,
Table[Sphere[i, 0.2],{i,1,Length[cData]}]
]
]

Again data set A looks quite random, while for both B and C we can see at certain angles something like ordered lattice planes.

Let' s get a list of distance from the origin and plot them

distancesA = Table[Norm[aData[[i]]],{i,1,Length[aData]}];distancesB = Table[Norm[bData[[i]]],{i,1,Length[bData]}];distancesC = Table[Norm[cData[[i]]],{i,1,Length[cData]}];

Short[distancesA]

Short[distancesB]

Short[distancesC]

{Histogram[distancesA,ImageSize->300],Histogram[distancesB,ImageSize->300],Histogram[distancesC,ImageSize->300]}

The distances seem to have a parabolic distribution when they are less than 5, as they must if the number of particles in any shell scales like r^2.

We need a decent measure of ρo, it doesn't really matter too much if we are off a little bit because it only will affect the scaling of g(r). Visually inspecting the data, the sphere of radius 5 is probably a good representative sampling range.

ρ0dataA = Length[Cases[distancesA,a_/;a≤5]]/(4/3 Pi 5^3);ρ0dataA//N

For data sset A, there are on average 2.2 particles in each 1 x 1 x 1 box.

ρ0dataB = Length[Cases[distancesB,a_/;a≤5]]/(4/3 Pi 5^3);ρ0dataB//N

For data sset B, there are on average ~2 particles in each 1 x 1 x 1 box.

ρ0dataC = Length[Cases[distancesC,a_/;a≤5]]/(4/3 Pi 5^3);ρ0dataC//N

For data set C, there is on average ~1 particle in each 1 x 1 x 1 box.

Sorting atoms according to their distances from zero, so that we can pick out the atoms around which we will do the g(r)

SortBy is a method to sort by an additional criterion

dataAOrdered = SortBy[aData,Norm];
dataAOrdered//Short
Last[dataAOrdered]

dataBOrdered = SortBy[bData,Norm];
dataBOrdered//Short
Last[dataBOrdered]

dataCOrdered = SortBy[cData,Norm];
dataCOrdered//Short
Last[dataCOrdered]

We will try to do the averaging for the first 10 atoms near the origin for each case, here are their positions:

nearestTenA = Take[dataAOrdered,10]

nearestTenB = Take[dataBOrdered,10]

nearestTenC = Take[dataCOrdered,10]

Now construct the ρ(r) function to count the number of atoms in each spherical shell around a given point

ρR[data_, origin_, rmax_, dr_,ρ0_] :=
Module[
{shiftedData, distances},
shiftedData = Table[data[[i]] -origin,{i,1,Length[data]}];
distances = Rest[Table[Norm[shiftedData[[i]]],{i,1,Length[data]}]];
Table[
{r + dr/2,Length[Cases[distances,a_/;(r ≤ a < r + dr)]]/(4 Pi (r + dr/2)^2 dr)/ρ0},
{r,0.0001,rmax,dr}
]
]

Plot  ρ(r) or g(r) function and try to make sense of the data sets

Plot g(r) using one point {0,0,0} as example for all three data sets, note that we can adjust delta_r for better visualization

ListLinePlot[ρR[aData,nearestTenA[[1]],5,.05,ρ0dataA],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Red,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set A: example g(r)"]

data set A looks quite random.

ListLinePlot[ρR[bData,nearestTenB[[1]],5,.02,ρ0dataB],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Red,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set B: example g(r)"]

data set B shows some periodicity, but not very clear what this might be.

ListLinePlot[ρR[cData,nearestTenC[[1]],5,.02,ρ0dataC],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Red,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set C: example g(r)"]

data set C appears to have peaks around 1 and Sqrt[2], which may be a clue that it might be the primitive cubic.

Let's work out a method for averaging over seveal atoms to make sure we are getting a representative g(r)

First we check and see if Mean will do what we want for several lists:

Now testing the Mean function

Mean[{{aa,bb, dd},{cc,ee,ff},{gg,ii,jj}}]

Obtain the averaged g(r) functions for each case, and plot them

severalρrA = Mean[
Table[ρR[aData,nearestTenA[[i]],5,.05,ρ0dataA],{i,1,10}]
];
Short[severalρrA]

severalρrB = Mean[
Table[ρR[bData,nearestTenB[[i]],5,.01,ρ0dataB],{i,1,10}]
];
Short[severalρrB]

severalρrC = Mean[
Table[ρR[cData,nearestTenC[[i]],5,.01,ρ0dataC],{i,1,10}]
];
Short[severalρrC]

dataAgRPlot = ListLinePlot[severalρrA,PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Blue,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set A: averaged g(r)"]

dataBgRPlot = ListLinePlot[severalρrB,PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Brown,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set B: averaged g(r)"]

dataCgRPlot = ListLinePlot[severalρrC,PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{ Red,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Data Set C: averaged g(r)"]

Generating perfect primitive cubic, body centered cubic (bcc) and face centered cubic (fcc) lattices

Generate cubic lattice

mySize=6;

Set size to be within [-6, 6].

cubicLattice =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];
cubicLattice//Short

Plot the primitive cubic lattice

plotCubic=ListPointPlot3D[cubicLattice,BoxRatios-> {1,1,1}]

Note that this cubic lattice has 1 atom within each 1 x 1 x 1 box.

Generate fcc lattice

Center 1 positions:

fccLatticeC1 =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0}+{1/2,1/2,0},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];
fccLatticeC1//Short

Center 2 positions:

fccLatticeC2 =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0}+{1/2,0,1/2},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];

Center 3 positions:

fccLatticeC3 =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0}+{0,1/2,1/2},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];

bccLatticeCtr =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0}+{1/2,1/2,1/2},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];

Obtain the fcc lattice and plot it.

fccLattice=Flatten[{fccLatticeC1,fccLatticeC2,fccLatticeC3,cubicLattice},1];

plotFCC=ListPointPlot3D[fccLattice,BoxRatios-> {1,1,1}]

Note that this fcc lattice has 4 atoms within each 1 x 1 x 1 box.

Generate bcc lattice

Center positions for bcc lattice:

bccLatticeCtr =
Flatten[
 Table[
i{0,0,1} + j {0,1,0} + k{1,0,0}+{1/2,1/2,1/2},
{i,-mySize,mySize},{j,-mySize,mySize},{k,-mySize,mySize}],
2
];

Obtain the bcc lattice and plot it.

bccLattice=Flatten[{bccLatticeCtr,cubicLattice},1];

plotBCC=ListPointPlot3D[bccLattice,BoxRatios-> {1,1,1}]

Note that this bcc lattice has 2 atoms within each 1 x 1 x 1 box.

Sorting to obtain the ordered data for perfect lattices with -6<x<6, -6<y<6, -6<z<6

dataCubicOrdered = SortBy[cubicLattice,Norm];
dataCubicOrdered//Short
Last[dataCubicOrdered]

dataFCCOrdered = SortBy[fccLattice,Norm];
dataFCCOrdered//Short
Last[dataFCCOrdered]

dataBCCOrdered = SortBy[bccLattice,Norm];
dataBCCOrdered//Short
Last[dataBCCOrdered]

Plotting the g(r) functions for the perfect lattices

cubicGrPlot = ListPlot[ρR[dataCubicOrdered,{0,0,0},5,.05,1],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{Green,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "Primitive Cubic Lattice g(r)"]

fccGrPlot = ListPlot[ρR[dataFCCOrdered,{0,0,0},5,.05,4],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{Pink,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "fcc Cubic Lattice g(r)"]

bccGrPlot = ListPlot[ρR[dataBCCOrdered,{0,0,0},5,.05,2],PlotRange-> All,AxesLabel-> {"Distance","g(r)"},ImageSize-> Large,PlotStyle->{Cyan,Thick},PlotMarkers-> Automatic,Filling-> Axis,BaseStyle-> Large,PlotLabel-> "bcc Cubic Lattice g(r)"]

Show[fccGrPlot,cubicGrPlot,bccGrPlot,PlotLabel-> "Comparison between cubic lattice and fcc lattice"]

We can clearly see that the signature peaks are different in each case!

Trying to infer possible lattice structure of the three data sets

Show[dataAgRPlot,cubicGrPlot,PlotLabel-> "Comparison between cubic lattice and data set A"]

data set A does not look like a primitive cubic lattice at all!

Show[dataAgRPlot,bccGrPlot,PlotLabel-> "Comparison between bcc lattice and data set A"]

data set A does not look like a bcc lattice at all!

Show[dataAgRPlot,fccGrPlot,PlotLabel-> "Comparison between fcc lattice and data set A"]

data set A does not look like a fcc lattice as well!

the g(r) goes to 1 for long range suggests that it is likely a amorphous material.

Show[dataBgRPlot,bccGrPlot,PlotLabel-> "Comparison between bcc lattice and data set B"]

data set B's "finger prints" seem to match well with a perfect bcc lattice... it is likely a thermalized bcc lattice!

Show[dataCgRPlot,cubicGrPlot,PlotLabel-> "Comparison between cubic lattice and data set C"]

data set C's "finger prints" seem to match well with a perfect primitive cubic lattice... it is probably a thermalized primitice cubic lattice!

Exploring with different plotting to visualize the data sets

One other check, let' s look at the data by slicing off a quadrant... data set C first

slicedcdata = Cases[cData,{x_,y_,z_}/;x>0 &&y>0 && z> 0];
Graphics3D[
GraphicsComplex[
slicedcdata,
Table[Sphere[i, 0.15],{i,1,Length[slicedcdata]}]
]
]

There is definitely some periodicity when we rotate the image a bit, but it is fairly difficult to tell what it is.  Other methods would include a fourier transform, which we will cover later...

Now trying to slice a cubic volume to see if it looks better...

slicedcdata = Cases[cData,{x_,y_,z_}/;x>0 && x< 4 &&y>0 && y< 4 && z> 0&& z<4];
Graphics3D[
GraphicsComplex[
slicedcdata,
Table[Sphere[i, 0.15],{i,1,Length[slicedcdata]}]
]
]

Rotate the image, the periodicity show up a lot better. If we know it is a cubic lattice, then it makes sense...

Try playing with data set B...

slicedcdata = Cases[bData,{x_,y_,z_}/;x>0 && x< 4 &&y>0 && y< 4 && z> 0&& z<4];
Graphics3D[
GraphicsComplex[
slicedcdata,
Table[Sphere[i, 0.3],{i,1,Length[slicedcdata]}]
]
]

Play with different sphere sizes and rotate, some angles we can see the close packed planes...

Try playing with data set A...

slicedcdata = Cases[aData,{x_,y_,z_}/;x>0 && x< 4 &&y>0 && y< 4 && z> 0&& z<4];
Graphics3D[
GraphicsComplex[
slicedcdata,
Table[Sphere[i, 0.3],{i,1,Length[slicedcdata]}]
]
]

Explore with different sphere sizes and rotate around, indeed it is quite random and hard to see any long range order.