1! 2! Copyright (C) 2001-2012 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!-------------------------------------------------------------- 9MODULE klist 10 ! 11 !! This module contains the variables related to the k-points. 12 ! 13 USE kinds, ONLY : DP 14 USE parameters, ONLY : npk 15 ! 16 IMPLICIT NONE 17 ! 18 SAVE 19 ! 20 CHARACTER (LEN=32) :: smearing 21 !! smearing type 22 REAL(DP) :: xk(3,npk) 23 !! coordinates of k points 24 REAL(DP) :: wk(npk) 25 !! weight of k points 26 REAL(DP) :: xqq(3) 27 !! coordinates of q point (used in the ACFDT part) 28 REAL(DP) :: degauss 29 !! smearing parameter 30 REAL(DP) :: nelec 31 !! number of electrons 32 REAL(DP) :: nelup=0.0_dp 33 !! number of spin-up electrons (if two_fermi_energies=t) 34 REAL(DP) :: neldw=0.0_dp 35 !! number of spin-dw electrons (if two_fermi_energies=t) 36 REAL(DP) :: tot_magnetization 37 !! nelup-neldw >= 0 (negative value means unspecified) 38 REAL(DP) :: tot_charge 39 !! total charge 40 REAL(DP) :: qnorm= 0.0_dp 41 !! |q|, used in phonon+US calculations only 42 INTEGER, ALLOCATABLE :: igk_k(:,:) 43 !! index of G corresponding to a given index of k+G 44 INTEGER, ALLOCATABLE :: ngk(:) 45 !! number of plane waves for each k point 46 ! 47 INTEGER :: nks 48 !! number of k points in this pool 49 INTEGER :: nkstot 50 !! total number of k points 51 INTEGER :: ngauss 52 !! type of smearing technique 53 LOGICAL :: lgauss 54 !! if .TRUE.: use gaussian broadening 55 LOGICAL :: ltetra 56 !! if .TRUE.: use tetrahedra 57 LOGICAL :: lxkcry=.FALSE. 58 !! if .TRUE.:k-pnts in cryst. basis accepted in input 59 LOGICAL :: two_fermi_energies 60 !! if .TRUE.: nelup and neldw set ef_up and ef_dw separately 61 ! 62CONTAINS 63 ! 64 !-------------------------------------------------------------- 65 SUBROUTINE init_igk( npwx, ngm, g, gcutw ) 66 !-------------------------------------------------------------- 67 !! Initialize indices \(\text{igk_k}\) and number of plane waves 68 !! per k-point: 69 ! 70 !! * \((k_{ik} + G)_i = k_{ik} + G_\text{igk}\); 71 !! * i = 1, \text{ngk}(\text{ik}); 72 !! * \text{igk} = \text{igk}_k(i,ik). 73 ! 74 INTEGER, INTENT (IN) :: npwx, ngm 75 REAL(DP), INTENT(IN) :: gcutw, g(3,ngm) 76 ! 77 REAL(DP), ALLOCATABLE :: gk(:) 78 INTEGER :: ik 79 ! 80 IF (.NOT.ALLOCATED(igk_k)) ALLOCATE( igk_k(npwx,nks) ) 81 IF (.NOT.ALLOCATED(ngk)) ALLOCATE( ngk(nks) ) 82 ! 83 ALLOCATE ( gk(npwx) ) 84 igk_k(:,:) = 0 85 ! 86 ! ... The following loop must NOT be called more than once in a run 87 ! ... or else there will be problems with variable-cell calculations 88 ! 89 DO ik = 1, nks 90 CALL gk_sort( xk(1,ik), ngm, g, gcutw, ngk(ik), igk_k(1,ik), gk ) 91 ENDDO 92 ! 93 DEALLOCATE( gk ) 94 ! 95 END SUBROUTINE init_igk 96 ! 97 ! 98 SUBROUTINE deallocate_igk( ) 99 ! 100 IF (ALLOCATED(ngk)) DEALLOCATE( ngk ) 101 IF (ALLOCATED(igk_k)) DEALLOCATE( igk_k ) 102 ! 103 END SUBROUTINE deallocate_igk 104 ! 105 ! 106END MODULE klist 107! 108! 109!-------------------------------------------------------------- 110MODULE lsda_mod 111 ! 112 !! It contains the variables needed for the LSDA calculation. 113 ! 114 USE kinds, ONLY : DP 115 USE parameters, ONLY : ntypx, npk 116 ! 117 IMPLICIT NONE 118 ! 119 SAVE 120 ! 121 LOGICAL :: lsda 122 !! true if lsda is active 123 REAL(DP) :: magtot 124 !! total magnetization 125 REAL(DP) :: absmag 126 !! total absolute magnetization 127 REAL(DP) :: starting_magnetization(ntypx) 128 !! the magnetization used to start with 129 INTEGER :: nspin 130 !! number of spin polarization: 2 if lsda, 1 other 131 INTEGER :: current_spin 132 !! spin of the current kpoint 133 INTEGER :: isk(npk) 134 !! for each k-point: 1=spin up, 2=spin down 135 ! 136END MODULE lsda_mod 137! 138! 139! 140MODULE rap_point_group 141 ! 142 USE kinds, ONLY : DP 143 ! 144 INTEGER :: code_group 145 !! The code of the point group 146 INTEGER :: nclass 147 !! The number of classes of the point group 148 INTEGER :: nelem(12) 149 !! The elements of each class 150 INTEGER :: elem(8,12) 151 !! Which elements in the smat list for each class 152 INTEGER :: which_irr(12) 153 !! For each class gives its position in the character table. 154 ! 155 COMPLEX(DP) :: char_mat(12,12) 156 !! the character tables: rap, class 157 ! 158 CHARACTER(LEN=15) :: name_rap(12) 159 !! the name of the representation 160 CHARACTER(LEN=3) :: ir_ram(12) 161 !! a string I, R or I+R for infrared, Raman, or 162 !! infrared+raman modes. 163 CHARACTER(LEN=11) :: gname 164 !! the name of the group 165 CHARACTER(LEN=5) :: name_class(12) 166 !! the name of the class 167 CHARACTER(LEN=55) :: elem_name(8,12)=' ' 168 !! the name of each symmetry in each class 169 ! 170END MODULE rap_point_group 171! 172! 173! 174MODULE rap_point_group_so 175 ! 176 USE kinds, ONLY : DP 177 ! 178 INTEGER :: nrap 179 !! The number of classes of the point group 180 INTEGER :: nelem_so(24) 181 !! The elements of each class 182 INTEGER :: elem_so(12,24) 183 !! Which elements in the smat list for each class 184 INTEGER :: has_e(12,24) 185 !! if -1 the smat is multiplied by -E 186 INTEGER :: which_irr_so(24) 187 !! For each class gives its position in the character table. 188 ! 189 COMPLEX(DP) :: char_mat_so(12,24) 190 !! the character tables 191 COMPLEX(DP) :: d_spin(2,2,48) 192 !! the rotation in spin space 193 ! 194 CHARACTER(len=15) :: name_rap_so(12) 195 !! the name of the representation 196 CHARACTER(len=5) :: name_class_so(24) 197 !! the name of the class 198 CHARACTER(len=5) :: name_class_so1(24) 199 !! the name of the class 200 CHARACTER(len=55) :: elem_name_so(12,24)=' ' 201 !! the name of each symmetry in each class 202 ! 203END MODULE rap_point_group_so 204! 205! 206! 207MODULE rap_point_group_is 208 ! 209 USE kinds, ONLY : DP 210 ! 211 INTEGER :: nsym_is 212 !! The number of operations of the invariant subgroup 213 INTEGER :: code_group_is 214 !! The code of the point invariant subgroup 215 ! 216 REAL(DP) :: ft_is(3,48) 217 !! The fractional transl. of the invariant subgroup 218 REAL(DP) :: sr_is(3,3,48) 219 !! The matrices of the invariant subgroup 220 ! 221 COMPLEX(DP) :: d_spin_is(2,2,48) 222 !! the rotation in spin space 223 ! 224 CHARACTER(LEN=45) :: sname_is(48) 225 !! name of the symmetries 226 CHARACTER(LEN=11) :: gname_is 227 !! the name of the invariant group 228 ! 229END MODULE rap_point_group_is 230! 231! 232! 233MODULE vlocal 234 ! 235 !! The variables needed for the local potential in reciprocal space. 236 ! 237 USE kinds, ONLY : DP 238 USE parameters, ONLY : ntypx 239 ! 240 SAVE 241 ! 242 COMPLEX(DP), ALLOCATABLE :: strf(:,:) 243 !! the structure factor 244 REAL(DP), ALLOCATABLE :: vloc(:,:) 245 !! the local potential for each atom type 246 REAL(DP) :: starting_charge(ntypx) 247 !! the atomic charge used to start with 248 ! 249END MODULE vlocal 250! 251! 252! 253MODULE wvfct 254 ! 255 !! The variables needed to compute the band structure 256 ! 257 USE kinds, ONLY : DP 258 ! 259 SAVE 260 ! 261 INTEGER :: npwx 262 !! maximum number of PW for wavefunctions 263 INTEGER :: nbndx 264 !! max number of bands use in iterative diag 265 INTEGER :: nbnd 266 !! number of bands 267 INTEGER :: npw 268 !! the number of plane waves 269 INTEGER :: current_k 270 !! the index of k-point under consideration 271 REAL(DP), ALLOCATABLE :: et(:,:) 272 !! eigenvalues of the hamiltonian 273 REAL(DP), ALLOCATABLE :: wg(:,:) 274 !! the weight of each k point and band 275 REAL(DP), ALLOCATABLE :: g2kin(:) 276 !! kinetic energy 277 INTEGER, ALLOCATABLE :: btype(:,:) 278 !! one if the corresponding state has to be 279 !! converged to full accuracy, zero otherwise 280 ! 281END MODULE wvfct 282! 283! 284! 285MODULE ener 286 ! 287 !! The variables needed to compute the energies 288 ! 289 USE kinds, ONLY : DP 290 ! 291 SAVE 292 ! 293 REAL(DP) :: etot 294 !! the total Kohn-Sham energy of the solid 295 REAL(DP) :: hwf_energy 296 !! this is the Harris-Weinert-Foulkes energy 297 REAL(DP) :: eband 298 !! the band energy 299 REAL(DP) :: deband 300 !! scf correction to have variational energy 301 REAL(DP) :: ehart 302 !! the Hartree energy 303 REAL(DP) :: etxc 304 !! the exchange and correlation energy 305 REAL(DP) :: vtxc 306 !! another exchange-correlation energy 307 REAL(DP) :: etxcc 308 !! the nlcc exchange and correlation 309 REAL(DP) :: ewld 310 !! the ewald energy 311 REAL(DP) :: elondon 312 !! the semi-empirical dispersion energy 313 REAL(DP) :: edftd3 314 !! the Grimme-d3 dispersion energy 315 REAL(DP) :: exdm 316 !! the XDM dispersion energy 317 REAL(DP) :: demet 318 !! variational correction ("-TS") for metals 319 REAL(DP) :: epaw 320 !! sum of one-center paw contributions 321 REAL(DP) :: ef 322 !! the Fermi energy 323 REAL(DP) :: ef_up 324 !! the Fermi energy up (if two_fermi_energies=.TRUE.) 325 REAL(DP) :: ef_dw 326 !! the Fermi energy down (if two_fermi_energies=.TRUE.) 327 ! 328END MODULE ener 329! 330! 331! 332MODULE force_mod 333 ! 334 !! The variables for the first derivative of the energy 335 ! 336 USE kinds, ONLY : DP 337 ! 338 SAVE 339 ! 340 REAL(DP), ALLOCATABLE :: force(:,:) 341 !! the force on each atom 342 REAL(DP) :: sumfor 343 !! norm of the gradient (forces) 344 REAL(DP) :: sigma(3,3) 345 !! the stress acting on the system 346 LOGICAL :: lforce 347 !! if .TRUE. compute the forces 348 LOGICAL :: lstres 349 !! if .TRUE. compute the stress 350 REAL(DP), ALLOCATABLE :: eigenval(:) 351 !! eigenvalues of the overlap matrix 352 COMPLEX(DP), ALLOCATABLE :: eigenvect(:,:) 353 !! eigenvectors of the overlap matrix 354 COMPLEX(DP), ALLOCATABLE :: overlap_inv(:,:) 355 !! overlap matrix (transposed): (O^{-1/2})^T 356 COMPLEX(DP), ALLOCATABLE :: doverlap_inv(:,:) 357 !! derivative of the overlap matrix (not transposed): d(O^{-1/2}) 358 COMPLEX (DP), ALLOCATABLE :: at_dy(:,:), at_dj(:,:) 359 !! derivatives of spherical harmonics and spherical Bessel functions (for atomic functions) 360 COMPLEX (DP), ALLOCATABLE :: us_dy(:,:), us_dj(:,:) 361 !! derivatives of spherical harmonics and spherical Bessel functions (for beta functions) 362 ! 363END MODULE force_mod 364! 365! 366! 367MODULE relax 368 ! 369 !! The variables used to control ionic relaxations 370 ! 371 USE kinds, ONLY : DP 372 ! 373 SAVE 374 ! 375 REAL(DP) :: epse = 0.0_dp 376 !! threshold on total energy 377 REAL(DP) :: epsf 378 !! threshold on forces 379 REAL(DP) :: epsp 380 !! threshold on pressure 381 REAL(DP) :: starting_scf_threshold 382 !! self-explanatory 383 ! 384END MODULE relax 385! 386! 387! 388MODULE cellmd 389 ! 390 !! The variables used to control cell relaxation. 391 ! 392 USE kinds, ONLY : DP 393 ! 394 SAVE 395 ! 396 REAL(DP) :: at_old(3,3) 397 !! the lattice vectors at the previous step 398 REAL(DP) :: omega_old 399 !! the cell volume at the previous step 400 REAL(DP) :: cell_factor=0.0_dp 401 !! maximum expected (linear) cell contraction 402 !! during relaxation/MD 403 INTEGER :: nzero 404 !! iteration # of last thermalization 405 INTEGER :: ntimes=-1 406 !! # of thermalization steps to be performed (-i=inf) 407 INTEGER :: ntcheck 408 !! # of steps between thermalizations 409 LOGICAL :: lmovecell 410 !! used in cell relaxation 411 ! 412 CHARACTER(LEN=2) :: calc=' ' 413 !! main switch for variable cell shape MD 414 !! see move_ions and vcsmd for allowed values 415 ! 416END MODULE cellmd 417! 418! 419! 420MODULE us 421 ! 422 !! These parameters are needed with the US pseudopotentials. 423 ! 424 USE kinds, ONLY : DP 425 ! 426 SAVE 427 ! 428 INTEGER :: nqxq 429 !! size of interpolation table 430 INTEGER :: nqx 431 !! number of interpolation points 432 REAL(DP), PARAMETER:: dq = 0.01D0 433 !! space between points in the pseudopotential tab. 434 REAL(DP), ALLOCATABLE :: qrad(:,:,:,:) 435 !! radial FT of Q functions 436 REAL(DP), ALLOCATABLE :: tab(:,:,:) 437 !! interpolation table for PPs 438 REAL(DP), ALLOCATABLE :: tab_at(:,:,:) 439 !! interpolation table for atomic wfc 440 LOGICAL :: spline_ps = .FALSE. 441 REAL(DP), ALLOCATABLE :: tab_d2y(:,:,:) 442 !! for cubic splines 443 ! 444END MODULE us 445! 446! 447! 448MODULE fixed_occ 449 ! 450 !! The quantities needed in calculations with fixed occupations. 451 ! 452 USE kinds, ONLY : DP 453 ! 454 SAVE 455 ! 456 REAL(DP), ALLOCATABLE :: f_inp(:,:) 457 !! the occupations for each spin 458 LOGICAL :: tfixed_occ 459 !! if .TRUE. the occupations are fixed. 460 LOGICAL :: one_atom_occupations 461 !! if .TRUE. the occupations are decided according to the projections of the 462 !! wavefunctions on the initial atomic wavefunctions (to be used only 463 !! for an isolated atom). 464 ! 465END MODULE fixed_occ 466! 467! 468! 469MODULE spin_orb 470 ! 471 !! Variables needed for calculations with spin-orbit 472 ! 473 USE kinds, ONLY : DP 474 USE upf_params, ONLY : lmaxx, lqmax 475 !! FIXME: rot_ylm could be dynamically allocated 476 ! 477 SAVE 478 ! 479 LOGICAL :: lspinorb 480 !! if .TRUE. this is a spin-orbit calculation 481 LOGICAL :: lforcet 482 !! if .TRUE. apply Force Theorem to calculate MAE 483 LOGICAL :: starting_spin_angle 484 !! if .TRUE. the initial wavefunctions are spin-angle functions. 485 LOGICAL :: domag 486 !! if .TRUE. magnetization is computed 487 COMPLEX (DP) :: rot_ylm(lqmax,lqmax) 488 !! transform real spherical harmonics into complex ones 489 COMPLEX (DP), ALLOCATABLE :: fcoef(:,:,:,:,:) 490 !! function needed to account for spinors. 491 ! 492END MODULE spin_orb 493! 494! 495! 496MODULE pwcom 497 ! 498 USE klist 499 USE lsda_mod 500 USE vlocal 501 USE wvfct 502 USE ener 503 USE force_mod 504 USE relax 505 USE cellmd 506 USE fixed_occ 507 USE spin_orb 508 ! 509END MODULE pwcom 510