File: vecfbz.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (55 lines) | stat: -rw-r--r-- 1,486 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55

! Copyright (C) 2008 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU Lesser General Public
! License. See the file COPYING for license details.

!BOP
! !ROUTINE: vecfbz
! !INTERFACE:
subroutine vecfbz(eps,bvec,vpl)
! !INPUT/OUTPUT PARAMETERS:
!   eps  : zero component tolerance (in,real)
!   bvec : reciprocal lattice vectors (in,real(3,3))
!   vpl  : input vector in lattice coordinates (inout,real(3))
! !DESCRIPTION:
!   Maps a vector in lattice coordinates to the first Brillouin zone. This is
!   done by first removing its integer components and then adding primitive
!   reciprocal lattice vectors until the shortest vector is found.
!
! !REVISION HISTORY:
!   Created September 2008 (JKD)
!EOP
!BOC
implicit none
! arguments
real(8), intent(in) :: eps,bvec(3,3)
real(8), intent(inout) :: vpl(3)
! local variables
integer i1,i2,i3,j1,j2,j3
real(8) v0(3),v1(3),v2(3),v3(3),t1,t2
! map vector to [0,1) interval
call r3frac(eps,vpl)
v0(:)=bvec(:,1)*vpl(1)+bvec(:,2)*vpl(2)+bvec(:,3)*vpl(3)
t1=v0(1)**2+v0(2)**2+v0(3)**2
j1=0; j2=0; j3=0
do i1=-1,0
  v1(:)=v0(:)+dble(i1)*bvec(:,1)
  do i2=-1,0
    v2(:)=v1(:)+dble(i2)*bvec(:,2)
    do i3=-1,0
      v3(:)=v2(:)+dble(i3)*bvec(:,3)
      t2=v3(1)**2+v3(2)**2+v3(3)**2
      if (t2.lt.t1+eps) then
        j1=i1; j2=i2; j3=i3
        t1=t2
      end if
    end do
  end do
end do
vpl(1)=vpl(1)+dble(j1)
vpl(2)=vpl(2)+dble(j2)
vpl(3)=vpl(3)+dble(j3)
return
end subroutine
!EOC