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