1! 2! Copyright (C) 2001-2007 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!---------------------------------------------------------------------------- 9SUBROUTINE deriv_drhoc( ngl, gl, omega, tpiba2, mesh, r, rab, rhoc, drhocg ) 10 !-------------------------------------------------------------------------- 11 !! Calculates the Fourier transform of \(d\text{Rho}_c/dG\). 12 ! 13 USE kinds 14 USE constants, ONLY : pi, fpi 15 ! 16 IMPLICIT NONE 17 ! 18 INTEGER :: ngl 19 !! input: the number of g shell 20 INTEGER :: mesh 21 !! input: the number of radial mesh points 22 REAL(DP), INTENT(IN) :: gl(ngl) 23 !! input: the number of G shells 24 REAL(DP), INTENT(IN) :: r(mesh) 25 !! input: the radial mesh 26 REAL(DP), INTENT(IN) :: rab(mesh) 27 !! input: the derivative of the radial mesh 28 REAL(DP), INTENT(IN) :: rhoc(mesh) 29 !! input: the radial core charge 30 REAL(DP), INTENT(IN) :: omega 31 !! input: the volume of the unit cell 32 REAL(DP), INTENT(IN) :: tpiba2 33 !! input: 2 times pi / alat 34 REAL(DP), INTENT(OUT) :: drhocg(ngl) 35 !! Fourier transform of d Rho_c/dG 36 ! 37 ! ... local variables 38 ! 39 REAL(DP) :: gx, rhocg1 40 ! the modulus of g for a given shell 41 ! the fourier transform 42 REAL(DP), ALLOCATABLE :: aux(:) 43 ! auxiliary memory for integration 44 INTEGER :: ir, igl, igl0 45 ! counter on radial mesh points 46 ! counter on g shells 47 ! lower limit for loop on ngl 48 ! 49 ! G=0 term 50 ! 51 IF (gl(1) < 1.0d-8) THEN 52 drhocg(1) = 0.0d0 53 igl0 = 2 54 ELSE 55 igl0 = 1 56 ENDIF 57 ! 58 ! G <> 0 term 59 ! 60!$omp parallel private(aux, gx, rhocg1) 61 ! 62 ALLOCATE( aux(mesh) ) 63!$omp do 64 DO igl = igl0, ngl 65 gx = SQRT( gl(igl) * tpiba2 ) 66 DO ir = 1, mesh 67 aux(ir) = r(ir)*rhoc(ir)*( r(ir) * COS(gx*r(ir)) / & 68 gx - SIN(gx*r(ir)) / gx**2 ) 69 ENDDO 70 CALL simpson( mesh, aux, rab, rhocg1 ) 71 drhocg(igl) = fpi / omega * rhocg1 72 ENDDO 73!$omp end do nowait 74 DEALLOCATE( aux ) 75 ! 76!$omp end parallel 77 ! 78 RETURN 79 ! 80END SUBROUTINE deriv_drhoc 81 82