File: allatoms.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 (62 lines) | stat: -rw-r--r-- 1,919 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

! 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: allatoms
! !INTERFACE:
subroutine allatoms
! !USES:
use modmain
use modxcifc
use modomp
! !DESCRIPTION:
!   Solves the Kohn-Sham-Dirac equations for each atom type in the solid and
!   finds the self-consistent radial wavefunctions, eigenvalues, charge
!   densities and potentials. The atomic densities can then be used to
!   initialise the crystal densities, and the atomic self-consistent potentials
!   can be appended to the muffin-tin potentials to solve for the core states.
!   Note that, irrespective of the value of {\tt xctype}, exchange-correlation
!   functional type 3 is used. See also {\tt atoms}, {\tt rhoinit},
!   {\tt gencore} and {\tt modxcifc}.
!
! !REVISION HISTORY:
!   Created September 2002 (JKD)
!   Modified for GGA, June 2007 (JKD)
!EOP
!BOC
implicit none
logical hybrid_
integer xcspin_,xcgrad_
integer is,nthd
real(8) hybridc_
character(512) xcdescr_
! allocatable arrays
real(8), allocatable :: rwf(:,:,:)
! allocate global species charge density and potential arrays
if (allocated(rhosp)) deallocate(rhosp)
allocate(rhosp(nrspmax,nspecies))
if (allocated(vrsp)) deallocate(vrsp)
allocate(vrsp(nrspmax,nspecies))
! get the exchange-correlation functional data
call getxcdata(xctsp,xcdescr_,xcspin_,xcgrad_,hybrid_,hybridc_)
call omp_hold(nspecies,nthd)
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(rwf) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
do is=1,nspecies
  allocate(rwf(nrspmax,2,nstspmax))
  call atom(solsc,ptnucl,spzn(is),nstsp(is),nsp(:,is),lsp(:,is),ksp(:,is), &
   occsp(:,is),xctsp,xcgrad_,nrsp(is),rsp(:,is),evalsp(:,is),rhosp(:,is), &
   vrsp(:,is),rwf)
  deallocate(rwf)
end do
!$OMP END DO
!$OMP END PARALLEL
call omp_free(nthd)
return
end subroutine
!EOC