Line data Source code
1 : MODULE m_grdrsvac
2 : use m_juDFT
3 : c----------------------------------------------------------------------
4 : c for a function defined on a vacuum layer
5 : c the in-plane derivatives are evaluated in real space
6 : c
7 : c based on 'rhzgrd' coded by t.asada. june,1995.
8 : c---------------------------------------------------------------------
9 : CONTAINS
10 0 : SUBROUTINE grdrsvac(
11 0 : > ro,bmat,xmax1,xmax2,ndvgrd,
12 0 : < dxro,dyro)
13 :
14 : c----------------------------------------------------------------------
15 : c
16 : c input:
17 : c ro(0:xmax1*xmax2-1)
18 : c any quantity stored in the 2-dim box (xmax1 x xmax2)
19 : c bmat
20 : c bravais matrix of reciprocal space
21 : c ndvgrd
22 : c number of points used when calculating derivative (3 <= ndvgrd <= 6)
23 : c
24 : c output:
25 : c dxro(0:xmax1*xmax2-1) , dxro(0:xmax1*xmax2-1)
26 : c (d ro / d x) , (d ro /d y) in non-internal coordinates
27 : c
28 : c----------------------------------------------------------------------
29 : c
30 : USE m_constants, ONLY : pimach
31 : IMPLICIT NONE
32 : c ..
33 : c .. Scalar arguments ..
34 : INTEGER, INTENT (IN) :: xmax1,xmax2,ndvgrd
35 : c ..
36 : c .. Array arguments ..
37 : REAL, INTENT (IN) :: ro(0:xmax1*xmax2-1)
38 : REAL, INTENT (IN) :: bmat(3,3)
39 : c ..
40 : c .. Array output ..
41 : REAL, INTENT (OUT) :: dxro(0:xmax1*xmax2-1)
42 : REAL, INTENT (OUT) :: dyro(0:xmax1*xmax2-1)
43 : c ..
44 : c .. Locals ..
45 : INTEGER :: xmax(2)
46 : INTEGER :: direction,xyz(2),x1,x2,i,ii(-3:2)
47 : REAL :: pi
48 : REAL :: dx
49 0 : REAL :: drointern(0:xmax1*xmax2-1,2)
50 :
51 0 : pi = pimach()
52 :
53 0 : xmax(1)= xmax1
54 0 : xmax(2)= xmax2
55 0 : DO i=1,2
56 0 : IF ( xmax(i) < 3 ) THEN
57 0 : CALL juDFT_error("grid to small",calledby="grdrsvac")
58 : END IF
59 : END DO
60 0 : IF ( (ndvgrd < 3) .or. (ndvgrd > 6) ) THEN
61 0 : CALL juDFT_error("ndvgrd notin [3,6]",calledby="grdrsvac")
62 : ENDIF
63 :
64 :
65 0 : DO direction=1,2
66 :
67 0 : dx= 1./REAL(xmax(direction))
68 :
69 0 : DO x1=0,xmax(1)-1
70 0 : DO x2=0,xmax(2)-1
71 :
72 0 : DO i= -3,2
73 0 : xyz(1)= x1
74 0 : xyz(2)= x2
75 0 : xyz(direction)= xyz(direction)+i
76 : ! make use of periodic boundary cond. in interstitial:
77 0 : IF ( xyz(direction) < 0 ) THEN
78 0 : xyz(direction)= xyz(direction)+xmax(direction)
79 : END IF
80 0 : IF ( xyz(direction) >= xmax(direction) ) THEN
81 0 : xyz(direction)= xyz(direction)-xmax(direction)
82 : END IF
83 : ! find coordinates in 1-dim array ro:
84 0 : ii(i)= xyz(2)*xmax(1) + xyz(1)
85 : END DO
86 :
87 0 : IF (ndvgrd.EQ.3) THEN
88 : drointern(ii(0),direction)=
89 : & df3( ro(ii(-1)),
90 0 : & ro(ii(0)),ro(ii(1)), dx)
91 0 : ELSEIF (ndvgrd.EQ.4) THEN
92 : drointern(ii(0),direction)=
93 : & df4( ro(ii(-1)),
94 0 : & ro(ii(0)),ro(ii(1)),ro(ii(2)), dx)
95 0 : ELSEIF (ndvgrd.EQ.5) THEN
96 : drointern(ii(0),direction)=
97 : & df5( ro(ii(-2)),ro(ii(-1)),
98 0 : & ro(ii(0)),ro(ii(1)),ro(ii(2)), dx)
99 0 : ELSEIF (ndvgrd.EQ.6) THEN
100 : drointern(ii(0),direction)=
101 : & df6( ro(ii(-3)),ro(ii(-2)),ro(ii(-1)),
102 0 : & ro(ii(0)),ro(ii(1)),ro(ii(2)), dx)
103 : ENDIF
104 :
105 : END DO
106 : END DO
107 :
108 : END DO
109 :
110 :
111 0 : DO i=0,xmax(1)*xmax(2)-1
112 :
113 : dxro(i)= bmat(1,1)*drointern(i,1)
114 0 : & + bmat(2,1)*drointern(i,2)
115 0 : dxro(i)= dxro(i)/(2.*pi)
116 : dyro(i)= bmat(1,2)*drointern(i,1)
117 0 : & + bmat(2,2)*drointern(i,2)
118 0 : dyro(i)= dyro(i)/(2.*pi)
119 :
120 : END DO
121 :
122 0 : END SUBROUTINE grdrsvac
123 : !--------------------------------------------------------------------
124 : ! Functions: formulae for 1st deriv.:
125 : !
126 0 : REAL FUNCTION df3(g1,f0,f1,d) ! three point formula
127 : REAL g1,f0,f1,d
128 0 : df3 = (-1*g1-0*f0+f1)/ (2*d)
129 0 : END FUNCTION df3
130 :
131 0 : REAL FUNCTION df4(g1,f0,f1,f2,d) ! four point formula
132 : REAL g1,f0,f1,f2,d
133 0 : df4 = (-2*g1-3*f0+6*f1-f2)/ (6*d)
134 0 : END FUNCTION df4
135 :
136 0 : REAL FUNCTION df5(g2,g1,f0,f1,f2,d) ! five point formula
137 : REAL g2,g1,f0,f1,f2,d
138 0 : df5 = (2*g2-16*g1-0*f0+16*f1-2*f2)/ (24*d)
139 0 : END FUNCTION df5
140 :
141 0 : REAL FUNCTION df6(g3,g2,g1,f0,f1,f2,d) ! six point formula
142 : REAL g3,g2,g1,f0,f1,f2,d
143 0 : df6 = (-4*g3+30*g2-120*g1+40*f0+60*f1-6*f2)/ (120*d)
144 0 : END FUNCTION df6
145 :
146 : !----------------------------------------------------------------------
147 : END MODULE m_grdrsvac
148 :
|
