File: gencfun.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 (68 lines) | stat: -rw-r--r-- 2,042 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
56
57
58
59
60
61
62
63
64
65
66
67
68

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

!BOP
! !ROUTINE: gencfun
! !INTERFACE:
subroutine gencfun
! !USES:
use modmain
! !DESCRIPTION:
!   Generates the smooth characteristic function. This is the function which is
!   0 within the muffin-tins and 1 in the intersitial region and is constructed
!   from radial step function form factors with $G<G_{\rm max}$. The form
!   factors are given by
!   $$ \tilde{\Theta}_i(G)=\begin{cases}
!    \frac{4\pi R_i^3}{3 \Omega} & G=0 \\
!    \frac{4\pi R_i^3}{\Omega}\frac{j_1(GR_i)}{GR_i} & 0<G\le G_{\rm max} \\
!    0 & G>G_{\rm max}\end{cases} $$
!   where $R_i$ is the muffin-tin radius of the $i$th species and $\Omega$ is
!   the unit cell volume. Therefore the characteristic function in $G$-space is
!   $$ \tilde{\Theta}({\bf G})=\delta_{G,0}-\sum_{ij}\exp(-i{\bf G}\cdot
!    {\bf r}_{ij})\tilde{\Theta}_i(G), $$
!   where ${\bf r}_{ij}$ is the position of the $j$th atom of the $i$th species.
!
! !REVISION HISTORY:
!   Created January 2003 (JKD)
!EOP
!BOC
implicit none
! local variables
integer is,ia,ig
real(8) v(3),t1
complex(8) z1
! allocatable arrays
complex(8), allocatable :: zfft(:)
! allocate global characteristic function arrays
if (allocated(cfunig)) deallocate(cfunig)
allocate(cfunig(ngtot))
if (allocated(cfunir)) deallocate(cfunir)
allocate(cfunir(ngtot))
cfunig(:)=0.d0
cfunig(1)=1.d0
! begin loop over species
do is=1,nspecies
! loop over atoms
  do ia=1,natoms(is)
    v(:)=atposc(:,ia,is)
    do ig=1,ngtot
! structure factor
      t1=vgc(1,ig)*v(1)+vgc(2,ig)*v(2)+vgc(3,ig)*v(3)
      z1=cmplx(cos(t1),-sin(t1),8)
! add to characteristic function in G-space
      cfunig(ig)=cfunig(ig)-ffacg(ig,is)*z1
    end do
  end do
end do
allocate(zfft(ngtot))
zfft(igfft(:))=cfunig(:)
! Fourier transform to real-space
call zfftifc(3,ngridg,1,zfft)
cfunir(:)=dble(zfft(:))
deallocate(zfft)
return
end subroutine
!EOC