Animation of the formation of the Si(111) 7x7 surface reconstruction
Matlab code
% Movie of the reconstruction of Si(111)7x7 starting with seven bulk atomic
% layers and allowing the top three to reconfigure to the 7x7
% reconstruction. The input data is generated by the program si7x7.m using
% atomic coordinates from the paper by K.D. Brommer et al., PRL, vol 68, pp
% 1355 - 1358 (1992).
clear; close all;
si7x7 = load('Si7x7coords.dat');
% Arrange data in descending order of their vertical position
si7x7 = sortrows(si7x7,3,'descend');
vid = VideoWriter('Si7x7.mp4','MPEG-4');
vid.Quality = 100;
vid.FrameRate = 30;
open(vid);
figure('units','pixels','position',[0 0 1920 1080],'ToolBar','none');
set(0,'defaultfigurecolor',[1 1 1]);
set(gca,'linewidth',7);
% Comment out one of the two below colormaps depending on your preference
% % Customized colormap from dark grey through red, orange, and
% % yellow in order to highlight different layer heights
% blackorange = customcolormap([0 0.35 0.88 1], ...
% {'#252525','#ff0000','#ff6a00','#ffff00'});
% blackorange = flip(blackorange);
% colormap(blackorange);
% Customized colormap from dark grey through the colors of the rainbow in
% order to highlight different layer heights
rainbow = customcolormap([0 1 2 3 4 5 6 7]/7, ...
{'#252525','#ff0000','#ff6a00','#ffff00','#00ff00','#00ffff','#0000ff','#7f00ff'});
rainbow = flip(rainbow);
colormap(rainbow);
colormap('jet');
c = jet;
rt3 = sqrt(3);
rt24 = sqrt(24);
[x,y,z] = sphere(10);
r = 0.25;
xmin = 0;
xmax = 14*rt3;
ymin = -7;
ymax = 7;
zmin = -3.5;
zmax = 0.5;
% Bulk atoms within the diamond reconstruction volume (7 x 7 x 3), composed
% of three layers of 7x7, plus six further bulk atomic layers that remain
% untouched. The total number of atoms is thus 7 x 7 x 9 = 441.
bulkCoords = zeros(size(441,1),4); % (x,y,z,index of associated (i.e., closest) surface atom)
indx = 1;
xshift = [0 0 -1 -1 -2 -2 -3 -3 -4]/rt3;
zshift = [1 4 5 8 9 12 13 16 17]/rt24;
for kk = 0:8 % z-layers
for jj = 0:6 % y-direction
for ii = 0:6 % x-direction
bulkCoords(indx,1) = ii*rt3/2 + jj*rt3/2 + xshift(kk+1);
bulkCoords(indx,2) = (jj - ii)/2;
bulkCoords(indx,3) = -zshift(kk+1);
indx = indx+1;
end
end
end
for ii = 1:147 % Go through 147 surface atoms
% iadx = difference between x-coordinate of iith surface atom and all
% x-coordinates of the bulk atoms. Same for y- and z-coordinates
iadx = (si7x7(ii,1) - bulkCoords(:,1));
iady = (si7x7(ii,2) - bulkCoords(:,2));
iadz = (si7x7(ii,3) - bulkCoords(:,3));
% Vector of the square of the interatomic distances between iith
% surface atom and all bulk atoms
iad = iadx.^2 + iady.^2 + iadz.^2;
[MMM,ind] = min(iad); % ind = index of bulk atom closest to surface atom ii
while (bulkCoords(ind,4) ~= 0) % This bulk atom has already been selected
iad(ind) = iad(ind)+100; % Make this iad so large it will be ignored in next search
[MMM,ind] = min(iad); % Repeat search procedure to find new minimum
end
bulkCoords(ind,4) = ii;
end
recMaxHt = max(si7x7(:,3)); % Maximum height of atoms in reconstruction
recMinHt = min(si7x7(:,3)); % Minimum height of atoms in reconstruction
for tt = 0:0.005:1 % Loop over 201 frames between bulk surface and reconstructed surface
newplot
tt % Monitor progress
tic % Monitor speed
hold off
for jj = -2:2 % multiple diamond blocks along y-axis
for kk = 0:2 % multiple diamond blocks along x-axis
recShiftx = -3.5*rt3*mod(jj+1,2)+kk*7*rt3;
recShifty = 3.5*jj;
hold on
for ii = 1:size(bulkCoords,1)
sindex = bulkCoords(ii,4); % Index of surface atom associated with iith bulk atom
% Determine the coordinates of atom as it moves from the
% bulk to the surface position
if (bulkCoords(ii,4) ~= 0)
xposNow = bulkCoords(ii,1) + tt*(si7x7(sindex,1) ...
- bulkCoords(ii,1)) + recShiftx;
yposNow = bulkCoords(ii,2) + tt*(si7x7(sindex,2) ...
- bulkCoords(ii,2)) + recShifty;
zposNow = bulkCoords(ii,3) + tt*(si7x7(sindex,3) ...
- bulkCoords(ii,3));
if (zposNow >= recMinHt)
% Index of colormap of surface atom at this point of
% the development of the reconstruction
colNowIndex = ...
1 + round(tt*255*(zposNow - recMinHt)/(recMaxHt - recMinHt));
colNow = c(colNowIndex,:);
else
colNowIndex = 1;
colNow = [25 25 25]/256;
end
else
xposNow = bulkCoords(ii,1) + recShiftx;
yposNow = bulkCoords(ii,2) + recShifty;
zposNow = bulkCoords(ii,3);
colNow = [25 25 25]/256;
end
% Initial positions of shown atoms avoids atoms popping in
% and out of video
xposIn = bulkCoords(ii,1) + recShiftx;
yposIn = bulkCoords(ii,2) + recShifty;
zposIn = bulkCoords(ii,3);
if (xposIn >= xmin-0.1) && (xposIn <= xmax+0.1) && ...
(yposIn >= ymin) && (yposIn <= ymax) && ...
(zposIn >= zmin) && (zposIn <= zmax)
surf(r*x+xposNow,r*y+yposNow,r*z+zposNow,...
'FaceAlpha',1,'FaceColor', ...
colNow,'LineStyle',...
'none','FaceLighting','gouraud',...
'DiffuseStrength',1,'AmbientStrength',0.56);
end
end
% Fade in remaining surface atoms not picked by bulk atoms
for ii = 148:151
xposNow = si7x7(ii,1) + recShiftx;
yposNow = si7x7(ii,2) + recShifty;
zposNow = si7x7(ii,3);
colNowIndex = ... % Index of colormap of surface atom
1 + round(255*(zposNow - recMinHt)/(recMaxHt - recMinHt));
colNow = c(colNowIndex,:);
if (xposNow >= xmin) && (xposNow <= xmax) && ...
(yposNow >= ymin) && (yposNow <= ymax) && ...
(zposNow >= zmin) && (zposNow <= zmax)
surf(r*x+xposNow,r*y+yposNow,r*z+zposNow,...
'FaceAlpha',tt,'FaceColor', ...
colNow,'LineStyle',...
'none','FaceLighting','gouraud',...
'DiffuseStrength',1,'AmbientStrength',0.56);
end
end
end
end
lp = [0.5 0 0.5];
light('Position',lp,'Style','infinite');
set(gca, 'Projection','perspective');
view(41,50)
axis off
axis equal
xlim([xmin-3.5*r xmax+3.5*r]);
ylim([ymin-3.2*r ymax+3.2*r]);
zlim([zmin-r zmax+r]);
frame = getframe(gcf);
writeVideo(vid,frame);
hold off
toc % Monitor speed
end
pauseMov(vid,77); % Pause for ca. 2.5 seconds to absorb the beauty
% Output the movie as an mpg file
close(vid)