Hex, Bugs and More Physics | Emre S. Tasci

I’ve written some functions to aid in checking the transformations as well as switching between the parameter ↔ vector representations.

latpar2vec([a,b,c,alpha,beta,gamma])
Given the a,b,c,α,β and γ, the function calculates the lattice vectors

function f=latpar2vec(X)
% function f=latpar2vec([a,b,c,alpha,beta,gamma])
% Calculates the lattice vectors from the given lattice paramaters.
  
a     = X(1);
b     = X(2);
c     = X(3);
alpha = X(4);
beta  = X(5);
gamma = X(6);
 
% Converting the angles to radians:
alpha = d2r(alpha);
beta  = d2r(beta);
gamma = d2r(gamma);
 
%Aligning a with the x axis:
ax = a;
ay = 0;
az = 0;
 
%Orienting b to lie on the x-y plane:
bz = 0;
 
% Now we have 6 unknowns for 6 equations..
bx = b*cos(gamma);
by = (b**2-bx**2)**.5;
cx = c*cos(beta);
cy = (b*c*cos(alpha)-bx*cx)/by;
cz = (c**2-cx**2-cy**2)**.5;
f=[ax ay az
   bx by bz
   cx cy cz];
endfunction

latvec2par([ax ay az; bx by bz; cx cy cz])
Calculates the lattice parameters a,b,c,α,β and γ from the lattice vectors

function f=latvec2par(x)
% function f=latvec2par(x)
%     ax ay az
% x = bx by bz
%     cx cy cz
%
% takes the unit cell vectors and calculates the lattice parameters
% Emre S. Tasci <e.tasci@tudelft.nl>   03/10/08
 
av = x(1,:);
bv = x(2,:);
cv = x(3,:);
 
a = norm(av);
b = norm(bv);
c = norm(cv);
 
alpha = acos((bv*cv')/(b*c))*180/pi;
beta  = acos((av*cv')/(a*c))*180/pi;
gamma = acos((av*bv')/(a*b))*180/pi;
 
f = [a b c alpha beta gamma];
endfunction

volcell([a,b,c,A,B,G])
Calculates the volume of the cell from the given lattice parameters, which is the determinant of the matrice build from the lattice vectors.

function f=volcell(X)
% function f=volcell([a,b,c,A,B,G])
% Calculate the cell volume from lattice parameters
%
% Emre S. Tasci, 09/2008
a=X(1);
b=X(2);
c=X(3);
A=X(4); % alpha
B=X(5); % beta
G=X(6); % gamma
 f=a*b*c*(1-cos(d2r(A))^2-cos(d2r(B))^2-cos(d2r(G))^2+2*cos(d2r(A))*cos(d2r(B))*cos(d2r(G)))^.5;
endfunction

Why’s there no “volcell” function for the unit cell vectors? You’re joking, right? (det(vector)) ! 🙂

Example

octave:13> % Define the unit cell for PtSn4 :
octave:13> A = latpar2vec([ 6.41900 11.35700  6.38800  90.0000  90.0000  90.0000 ])
A =
 
    6.41900    0.00000    0.00000
    0.00000   11.35700    0.00000
    0.00000    0.00000    6.38800
 
octave:14> % Cell volume :
octave:14> Apar = [ 6.41900 11.35700  6.38800  90.0000  90.0000  90.0000 ]
Apar =
 
    6.4190   11.3570    6.3880   90.0000   90.0000   90.0000
 
octave:15> % Define the unit cell for PtSn4 :
octave:15> Apar=[ 6.41900 11.35700  6.38800  90.0000  90.0000  90.0000 ]
Apar =
 
    6.4190   11.3570    6.3880   90.0000   90.0000   90.0000
 
octave:16> % Cell volume :
octave:16> Avol = volcell (Apar)
Avol =  465.69
 
octave:17> % Calculate the lattice vectors :
octave:17> A = latpar2vec (Apar)
A =
 
    6.41900    0.00000    0.00000
    0.00000   11.35700    0.00000
    0.00000    0.00000    6.38800
 
octave:18> % Verify the volume :
octave:18> det(A)
ans =  465.69
 
octave:19> % Define the transformation matrix :
octave:19> R = [ 0 0 -1 ; -1 0 0 ; .5 .5 0]
R =
 
   0.00000   0.00000  -1.00000
  -1.00000   0.00000   0.00000
   0.50000   0.50000   0.00000
 
octave:21> % The reduced unit cell volume will be half of the original one as is evident from :
octave:21> det(R)
ans =  0.50000
 
octave:22> % Do the transformation :
octave:22> N = R*A
N =
 
  -0.00000  -0.00000  -6.38800
  -6.41900   0.00000   0.00000
   3.20950   5.67850   0.00000
 
octave:23> % The reduced cell parameters :
octave:23> Npar = latvec2par (N)
Npar =
 
     6.3880     6.4190     6.5227   119.4752    90.0000    90.0000
 
octave:24> % And the volume :
octave:24> det(N), volcell (Npar)
ans =  232.84
ans =  232.84