1! 2! Copyright (C) 2001 PWSCF 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!---------------------------------------------------------------------------- 9! TB 10! setup of the gate, search for 'TB' 11!---------------------------------------------------------------------------- 12! 13!---------------------------------------------------------------------- 14SUBROUTINE setlocal 15 !---------------------------------------------------------------------- 16 !! This routine computes the local potential in real space vltot(ir). 17 ! 18 USE io_global, ONLY : stdout 19 USE kinds, ONLY : DP 20 USE constants, ONLY : eps8 21 USE ions_base, ONLY : zv, ntyp => nsp 22 USE cell_base, ONLY : omega 23 USE extfield, ONLY : tefield, dipfield, etotefield, gate, & 24 etotgatefield !TB 25 USE gvect, ONLY : igtongl, gg 26 USE scf, ONLY : rho, v_of_0, vltot 27 USE vlocal, ONLY : strf, vloc 28 USE fft_base, ONLY : dfftp 29 USE fft_interfaces, ONLY : invfft 30 USE gvect, ONLY : ngm 31 USE control_flags, ONLY : gamma_only 32 USE mp_bands, ONLY : intra_bgrp_comm 33 USE mp, ONLY : mp_sum 34 USE martyna_tuckerman, ONLY : wg_corr_loc, do_comp_mt 35 USE esm, ONLY : esm_local, esm_bc, do_comp_esm 36 USE qmmm, ONLY : qmmm_add_esf 37 USE Coul_cut_2D, ONLY : do_cutoff_2D, cutoff_local 38 ! 39 IMPLICIT NONE 40 ! 41 COMPLEX(DP), ALLOCATABLE :: aux(:), v_corr(:) 42 ! auxiliary variable 43 INTEGER :: nt, ng 44 ! counter on atom types 45 ! counter on g vectors 46 ! 47 ALLOCATE( aux(dfftp%nnr) ) 48 aux(:) = (0.d0,0.d0) 49 ! 50 IF (do_comp_mt) THEN 51 ALLOCATE( v_corr(ngm) ) 52 CALL wg_corr_loc( omega, ntyp, ngm, zv, strf, v_corr ) 53 aux(dfftp%nl(:)) = v_corr(:) 54 DEALLOCATE( v_corr ) 55 ENDIF 56 ! 57 DO nt = 1, ntyp 58 DO ng = 1, ngm 59 aux(dfftp%nl(ng)) = aux(dfftp%nl(ng)) + vloc(igtongl(ng),nt) & 60 * strf(ng,nt) 61 ENDDO 62 ENDDO 63 ! 64 IF (gamma_only) THEN 65 DO ng = 1, ngm 66 aux(dfftp%nlm(ng)) = CONJG( aux(dfftp%nl(ng)) ) 67 ENDDO 68 ENDIF 69 ! 70 IF ( do_comp_esm .AND. ( esm_bc .NE. 'pbc' ) ) THEN 71 ! 72 ! ... Perform ESM correction to local potential 73 ! 74 CALL esm_local( aux ) 75 ! 76 ENDIF 77 ! 78 ! 2D: re-add the erf/r function 79 IF ( do_cutoff_2D ) THEN 80 ! 81 ! ... re-add the CUTOFF fourier transform of erf function 82 ! 83 CALL cutoff_local( aux ) 84 ! 85 ENDIF 86 ! 87 ! ... v_of_0 is (Vloc)(G=0) 88 ! 89 v_of_0 = 0.0_DP 90 IF (gg(1) < eps8) v_of_0 = DBLE( aux(dfftp%nl(1)) ) 91 ! 92 CALL mp_sum( v_of_0, intra_bgrp_comm ) 93 ! 94 ! ... aux = potential in G-space . FFT to real space 95 ! 96 CALL invfft( 'Rho', aux, dfftp ) 97 ! 98 vltot(:) = DBLE( aux(:) ) 99 ! 100 ! ... If required add an electric field to the local potential 101 ! 102 IF ( tefield .AND. ( .NOT. dipfield ) ) & 103 CALL add_efield( vltot, etotefield, rho%of_r, .TRUE. ) 104 ! 105 ! TB 106 ! if charged plate, call add_gatefield and add the linear potential, 107 ! together with the background charge 108 IF (gate) CALL add_gatefield( vltot, etotgatefield, .TRUE., .TRUE. ) 109 ! 110 ! ... Add the electrostatic field generated by MM atoms 111 ! in a QM/MM calculation to the local potential 112 ! 113 CALL qmmm_add_esf( vltot, dfftp ) 114 ! 115 ! ... Save vltot for possible modifications in plugins 116 ! 117 CALL plugin_init_potential( vltot ) 118 ! 119 DEALLOCATE( aux ) 120 ! 121 ! 122 RETURN 123 ! 124END SUBROUTINE setlocal 125 126