Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_wann_plot_vac
8 : use m_juDFT
9 : c**************************************************************
10 : c Calculates the lattice periodic part of the Bloch function
11 : c in the vacuum. Used for plotting.
12 : c FF, Sep. '06
13 : c
14 : c***************************************************************
15 : CONTAINS
16 0 : SUBROUTINE wann_plot_vac(point,
17 : > z1,nmzd,nv2d,n3d,nvac,
18 : > nmz,delz,bmat,bbmat,
19 : > evac,bkpt,vz,
20 : > jspin,
21 0 : > k1,k2,k3,nvd,
22 0 : > nbasfcn,neigd,nv,omtil,nslibd,ac,bc,u,ue,
23 : < value)
24 :
25 : use m_constants, only : pimach
26 : implicit none
27 :
28 : c .. scalar Arguments..
29 : integer, intent (in) :: nmzd,nv2d,n3d,nslibd
30 : integer, intent (in) :: nmz,nvac
31 : integer, intent (in) :: jspin,nvd
32 : integer, intent (in) :: nbasfcn,neigd
33 : real, intent (in) :: delz,z1,omtil,point(3)
34 :
35 : c ..array arguments..
36 : real, intent (in) :: bkpt(3),evac(2)
37 : integer, intent (in) :: nv
38 : real, intent (in) :: vz(nmzd,2),bmat(3,3),bbmat(3,3)
39 : integer, intent (in) :: k1(nvd),k2(nvd),k3(nvd)
40 : complex, intent (out) :: value
41 :
42 : c ..basis wavefunctions in the vacuum
43 : complex, intent(in) :: ac(nv2d),bc(nv2d)
44 : real, intent(in) :: u(nmzd,nv2d),ue(nmzd,nv2d)
45 :
46 :
47 : c ..local scalars..
48 : real wronk,arg,zks,tpi,vz0(2),scale,evacp,ev,const
49 : real uu,ud,du,dd,xx(nmz)
50 : integer i,m,l,j,k,n,nv2,nv2_b,ivac,n2,n2_b,sign,ik
51 : integer lprime,np1
52 : complex av,bv,ic,c_1
53 0 : integer, allocatable :: kvac1(:),kvac2(:),map2(:)
54 : complex value1
55 : c ..intrinsic functions..
56 : intrinsic aimag,cmplx,conjg,real,sqrt
57 :
58 0 : call timestart("wann_plot_vac")
59 0 : allocate (kvac1(nv2d),kvac2(nv2d),map2(nvd))
60 :
61 0 : tpi = 2 * pimach() ; ic = cmplx(0.,1.)
62 :
63 0 : np1 = nmz + 1
64 :
65 : c.. determining the indexing array (in-plane stars)
66 : c.. for the k-point
67 :
68 0 : wronk = 2.0
69 : const = 1.0 / ( sqrt(omtil)*wronk )
70 :
71 0 : n2 = 0
72 :
73 0 : do 40 k = 1,nv
74 0 : do 30 j = 1,n2
75 0 : if ( k1(k).eq.kvac1(j) .and.
76 : + k2(k).eq.kvac2(j) ) then
77 0 : map2(k) = j
78 0 : goto 40
79 : endif
80 0 : 30 continue
81 0 : n2 = n2 + 1
82 0 : if (n2>nv2d) CALL juDFT_error("wann_plot_vac: map",calledby
83 0 : + ="wann_plot_vac")
84 0 : kvac1(n2) = k1(k)
85 0 : kvac2(n2) = k2(k)
86 0 : map2(k) = n2
87 0 : 40 continue
88 :
89 :
90 0 : nv2 = n2
91 :
92 : c.. the body of the routine
93 :
94 :
95 :
96 0 : value=cmplx(0.0,0.0)
97 0 : value1=cmplx(0.0,0.0)
98 : c print*,"difference=",(abs(point(3))-z1)/delz
99 0 : i=(abs(point(3))-z1)/delz +1
100 :
101 :
102 :
103 0 : if (i.gt.nmz) then
104 0 : i=nmz
105 0 : print*,"i.gt.nmz in wann_plot_vac"
106 : endif
107 :
108 0 : do l = 1,nv2 !calculation for i
109 :
110 :
111 : arg=(kvac1(l)*bmat(1,1)+kvac2(l)*bmat(2,1))*point(1)+
112 0 : + (kvac1(l)*bmat(1,2)+kvac2(l)*bmat(2,2))*point(2)
113 0 : c_1=cmplx(cos(arg),sin(arg))
114 0 : value =value+ (u(i,l)*ac(l)+ue(i,l)*bc(l))*c_1
115 : c print*,"value=",value
116 0 : if (real(value).gt.10.or.real(value).lt.-10)then
117 0 : print*,"value=",value
118 0 : print*,"i=",i
119 0 : print*,"u(i,l)=",u(i,l)
120 0 : print*,"ac(l)=",ac(l)
121 0 : print*,"bc(l)=",bc(l)
122 0 : print*,"ue(i,l)=",ue(i,l)
123 : endif
124 : enddo ! l
125 :
126 :
127 0 : i=i+1
128 0 : do l = 1,nv2
129 :
130 : arg=(kvac1(l)*bmat(1,1)+kvac2(l)*bmat(2,1))*point(1)+
131 0 : + (kvac1(l)*bmat(1,2)+kvac2(l)*bmat(2,2))*point(2)
132 0 : c_1=cmplx(cos(arg),sin(arg))
133 0 : value1 =value1+ (u(i,l)*ac(l)+ue(i,l)*bc(l))*c_1
134 :
135 : enddo ! l
136 :
137 0 : value=(value1-value)*((abs(point(3))-z1)/delz+2-i)+value
138 :
139 0 : deallocate (kvac1,kvac2,map2 )
140 0 : call timestop("wann_plot_vac")
141 0 : END SUBROUTINE wann_plot_vac
142 : END MODULE m_wann_plot_vac
|