[Pw_forum] Coordinate transorms: Crystal to Cartesian and back

[Kostyantyn Borysenko]

Hi everyone,
I'm having problems with the math behind some subroutines. I am talking
about PH/trntnsc.f90 and PW/cryst_to_car.f90. 
These subroutines are supposed to convert vectors (cryst_to_car.f90) or
matrices (trntnsc.f90) from crystal coordinates to Cartesian (flag=1) or
vice versa (flag=-1). And it seems like a very straightforward task - just a
little bit of vector algebra. But when I derive these transformations myself
they don't check out. Here is what I've got.
 
By definition of reciprocal vectors we have the following property for
dot-product:
(a(i),b(j))= 2*pi*delta(i,j),                (1a)
where a(i) - lattice vectors, b(j) - reciprocal lattice vectors. In QE these
vectors are stored in matrices at(i,j) and bg(i,j) respectively. Matrix
elements of at are in units a_0 and matrix elements of bg are in units
(2*pi/a_0). So the condition (1a) can be rewritten as follows:
 SUM(k){at(k,i)*bg(k,j)} = delta(i,j).                (1b)
In matrix form it can be represented as
B' * A = I, or B' = inv(A),                    (2)
where B' is a transposed matrix B and inv(A) is an inversed matrix A, and I
is the identity matrix. I checked these matrices for the case of graphene
and property (2) indeed holds true.
 
Next, it's easy to show that transformation of an arbitrary vector v from
crystal to Cartesian coordinates is described by the matrix equation:
v_cart = A*v_cryst,                (3a)
where A is a matrix at(i,j) mentioned above, which contains lattice vector
components. If you multiply equation (3a) from the left by inv(A) and use
the property inv(A)=B' (see eq. (2)) you will get the reverse
transformation:
v_cryst = B' * v_cart                (4a)
 
Obviously, for an arbitrary matrix W the transformation equations will be
W_cart = A' *W_cryst * A                (3b)
W_cryst = B * W_cart * B'                (4b)
 
But that's not what I see in the code. In case of subroutine trntnsc.f90
matrices A and B switch places in (3b) and (4b). And in cryst_to_car.f90
everything seems fine but only if you use at(i,j) as a transformation matrix
when flag=1 and bg(i,j) if flag=-1. However, when it's called in subroutine
lint.f90 it looks like this:
 
call cryst_to_cart (nkh,xp,at,-1)
 
Like I said, I think the correct result can be obtained only in two cases:
call lint(nkh,xp,at,+1)   - Crystal -> Cartesian
call lint(nkh,xp,bg,-1)   - Cartesian -> Crystal 
 
Maybe I missed something but this seems like a pretty simple problem from
vector algebra. Any ideas?
 
 
Best regards,
 
Kostyantyn Borysenko
Department of Electrical and Computer Engineering 



North Carolina State University


-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://www.democritos.it/pipermail/pw_forum/attachments/20090522/d4d74362/attachment-0001.htm