[Pw_forum] Coordinate transorms: Crystal to Cartesian and back
Kostyantyn Borysenko
kboryse at ncsu.edu
Fri May 22 07:26:54 CEST 2009
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://lists.quantum-espresso.org/pipermail/users/attachments/20090522/d4d74362/attachment.html>
More information about the users
mailing list