1! 2! Copyright (C) 2002-2020 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!--------------------------------------------- 9! TB 10! included gate related stuff, search 'TB' 11!--------------------------------------------- 12! 13!=----------------------------------------------------------------------------=! 14! 15MODULE input_parameters 16! 17!=----------------------------------------------------------------------------=! 18! 19! this module contains 20! 1) the definitions of all input parameters 21! (both those read from namelists and those read from cards) 22! 2) the definitions of all namelists 23! 3) routines that allocate data needed in input 24! Note that all values are initialized, but the default values should be 25! set in the appropriate routines contained in module "read_namelists" 26! The documentation of input variables can be found in Doc/INPUT_PW.* 27! (for pw.x) or in Doc/INPUT_CP (for cp.x) 28! Originally written by Carlo Cavazzoni for FPMD 29! 30!=----------------------------------------------------------------------------=! 31 ! 32 USE kinds, ONLY : DP 33 USE parameters, ONLY : nsx, natx, sc_size 34 USE wannier_new,ONLY : wannier_data 35 ! 36 IMPLICIT NONE 37 ! 38 SAVE 39 ! 40!=----------------------------------------------------------------------------=! 41! BEGIN manual 42! 43! 44! * DESCRIPTION OF THE INPUT FILE 45! (to be given as standard input) 46! 47! The input file has the following layout: 48! 49! &CONTROL 50! control_parameter_1, 51! control_parameter_2, 52! ....... 53! control_parameter_Lastone 54! / 55! &SYSTEM 56! sistem_parameter_1, 57! sistem_parameter_2, 58! ....... 59! sistem_parameter_Lastone 60! / 61! &ELECTRONS 62! electrons_parameter_1, 63! electrons_parameter_2, 64! ....... 65! electrons_parameter_Lastone 66! / 67! &IONS 68! ions_parameter_1, 69! ions_parameter_2, 70! ....... 71! ions_parameter_Lastone 72! / 73! &CELL 74! cell_parameter_1, 75! cell_parameter_2, 76! ....... 77! cell_parameter_Lastone 78! / 79! ATOMIC_SPECIES 80! slabel_1 mass_1 pseudo_file_1 81! slabel_2 mass_2 pseudo_file_2 82! ..... 83! ATOMIC_POSITIONS 84! alabel_1 px_1 py_1 pz_1 85! alabel_2 px_2 py_2 pz_2 86! ..... 87! CARD_3 88! .... 89! CARD_N 90! 91! -- end of input file -- 92! 93!=----------------------------------------------------------------------------=! 94! CONTROL Namelist Input Parameters 95!=----------------------------------------------------------------------------=! 96! 97 CHARACTER(len=80) :: title = ' ' 98 ! a string describing the current job 99 100 CHARACTER(len=80) :: calculation = 'none' 101 ! Specify the type of the simulation 102 ! See below for allowed values 103 CHARACTER(len=80) :: calculation_allowed(15) 104 DATA calculation_allowed / 'scf', 'nscf', 'relax', 'md', 'cp', & 105 'vc-relax', 'vc-md', 'vc-cp', 'bands', 'neb', 'smd', 'cp-wf', & 106 'vc-cp-wf', 'cp-wf-nscf', 'ensemble'/ 107 CHARACTER(len=80) :: verbosity = 'default' 108 ! define the verbosity of the code output 109 CHARACTER(len=80) :: verbosity_allowed(6) 110 DATA verbosity_allowed / 'debug', 'high', 'medium', 'default', & 111 'low', 'minimal' / 112 113 CHARACTER(len=80) :: restart_mode = 'restart' 114 ! specify how to start/restart the simulation 115 CHARACTER(len=80) :: restart_mode_allowed(3) 116 DATA restart_mode_allowed / 'from_scratch', 'restart', 'reset_counters' / 117 118 INTEGER :: nstep = 10 119 ! number of simulation steps, see "restart_mode" 120 121 INTEGER :: iprint = 10 122 ! number of steps/scf iterations between successive writings 123 ! of relevant physical quantities to standard output 124 INTEGER :: max_xml_steps = 0 125 ! max number of steps between successive appending of an xml step 126 ! in the xml data file, default 0 means all steps are printed. 127 INTEGER :: isave = 100 128 ! number of steps between successive savings of 129 ! information needed to restart the run (see "ndr", "ndw") 130 ! used only in CP 131 132 LOGICAL :: tstress = .true. 133 ! .TRUE. calculate the stress tensor 134 ! .FALSE. do not calculate the stress tensor 135 136 LOGICAL :: tprnfor = .true. 137 ! .TRUE. calculate the atomic forces 138 ! .FALSE. do not calculate the atomic forces 139 140 REAL(DP) :: dt = 1.0_DP 141 ! time step for molecular dynamics simulation, in atomic units 142 ! CP: 1 a.u. of time = 2.4189 * 10^-17 s, PW: twice that much 143 ! Typical values for CP simulations are between 1 and 10 a.u. 144 ! For Born-Oppenheimer simulations, larger values can be used, 145 ! since it mostly depends only upon the mass of ions. 146 147 INTEGER :: ndr = 50 148 ! Fortran unit from which the code reads the restart file 149 150 INTEGER :: ndw = 50 151 ! Fortran unit to which the code writes the restart file 152 153 CHARACTER(len=256) :: outdir = './' 154 ! specify the directory where the code opens output and restart 155 ! files. When possible put this directory in the fastest available 156 ! filesystem ( not NFS! ) 157 158 CHARACTER(len=256) :: prefix = 'prefix' 159 ! specify the prefix for the output file, if not specified the 160 ! files are opened as standard fortran units. 161 162 CHARACTER(len=256) :: pseudo_dir = './' 163 ! specify the directory containing the pseudopotentials 164 165 REAL(DP) :: refg = 0.05_DP 166 ! Accurancy of the interpolation table, interval between 167 ! table values in Rydberg 168 169 CHARACTER(len=256) :: wfcdir = 'undefined' 170 ! scratch directory that is hopefully local to the node 171 ! to store large, usually temporary files. 172 173 REAL(DP) :: max_seconds = 1.0E+7_DP 174 ! smoothly terminate program after the specified number of seconds 175 ! this parameter is typically used to prevent an hard kill from 176 ! the queuing system. 177 178 REAL(DP) :: ekin_conv_thr = 1.0E-5_DP 179 ! convergence criterion for electron minimization 180 ! this criterion is met when "ekin < ekin_conv_thr" 181 ! convergence is achieved when all criteria are met 182 183 REAL(DP) :: etot_conv_thr = 1.0E-4_DP 184 ! convergence criterion for ion minimization 185 ! this criterion is met when "etot(n+1)-etot(n) < etot_conv_thr", 186 ! where "n" is the step index, "etot" the DFT energy 187 ! convergence is achieved when all criteria are met 188 189 REAL(DP) :: forc_conv_thr = 1.0E-3_DP 190 ! convergence criterion for ion minimization 191 ! this criterion is met when "MAXVAL(fion) < forc_conv_thr", 192 ! where fion are the atomic forces 193 ! convergence is achieved when all criteria are met 194 195 CHARACTER(len=80) :: disk_io = 'default' 196 ! Specify the amount of I/O activities 197 198 LOGICAL :: tefield = .false. 199 ! if .TRUE. a sawtooth potential simulating a finite electric field 200 ! is added to the local potential = only used in PW 201 202 ! TB - added gate also below to the namelist 203 LOGICAL :: gate = .FALSE. 204 ! if .TRUE. a charged plate in charged systems is added with a 205 ! total charge which is opposite to the charge of the system 206 207 LOGICAL :: tefield2 = .false. 208 ! if .TRUE. a second finite electric field is added to the local potential 209 ! only used in CP 210 211 LOGICAL :: lelfield = .false. 212 ! if .TRUE. a static homogeneous electric field is present 213 ! via the modern theory of polarizability - differs from tefield! 214 215 LOGICAL :: lorbm = .false. 216 ! if .TRUE. an orbital magnetization is computed (Kubo terms) 217 218 LOGICAL :: dipfield = .false. 219 ! if .TRUE. the dipole field is subtracted 220 ! only used in PW for surface calculations 221 222 LOGICAL :: lberry = .false. 223 ! if .TRUE., use modern theory of the polarization 224 225 INTEGER :: gdir = 0 226 ! G-vector for polarization calculation ( related to lberry ) 227 ! only used in PW 228 229 INTEGER :: nppstr = 0 230 ! number of k-points (parallel vector) ( related to lberry ) 231 ! only used in PW 232 233 INTEGER :: nberrycyc = 1 234 !number of covergence cycles on electric field 235 236 LOGICAL :: wf_collect = .true. 237 ! This flag controls the way wavefunctions are stored to disk: 238 ! .TRUE. collect wavefunctions from all processors, store them 239 ! into a single restart file on a single processors 240 ! .FALSE. do not collect wavefunctions, store them into distributed 241 ! files - NO LONGER IMPLEMENTED SINCE v.6.3 242 243 LOGICAL :: saverho = .true. 244 ! This flag controls the saving of charge density in CP codes: 245 ! .TRUE. save charge density to restart dir 246 ! .FALSE. do not save charge density 247 248 LOGICAL :: tabps = .false. ! for ab-initio pressure and/or surface 249 ! calculations 250 251 LOGICAL :: lkpoint_dir = .false. ! obsolete, for compatibility 252 253 LOGICAL :: use_wannier = .false. ! use or not Wannier functions 254 255 LOGICAL :: lecrpa = .FALSE. 256 ! if true symmetry in scf run is neglected for RPA Ec calculation 257 ! 258 259 LOGICAL :: tqmmm = .FALSE. ! QM/MM coupling. enabled if .true. 260 261 CHARACTER(len=256) :: vdw_table_name = ' ' 262 263 CHARACTER(len=10) :: point_label_type='SC' 264 265 CHARACTER(len=80) :: memory = 'default' 266 ! controls memory usage 267 CHARACTER(len=80) :: memory_allowed(3) 268 DATA memory_allowed / 'small', 'default', 'large' / 269 ! if memory = 'small' then QE tries to use (when implemented) algorithms using less memory, 270 ! even if they are slower than the default 271 ! if memory = 'large' then QE tries to use (when implemented) algorithms using more memory 272 ! to enhance performance. 273 274 LOGICAL :: lfcpopt = .FALSE. ! FCP optimisation switch 275 LOGICAL :: lfcpdyn = .FALSE. ! FCP thermostat enabled if .true. 276 ! 277 ! location of xml input according to xsd schema 278 CHARACTER(len=256) :: input_xml_schema_file = ' ' 279 280 NAMELIST / control / title, calculation, verbosity, restart_mode, & 281 nstep, iprint, isave, tstress, tprnfor, dt, ndr, ndw, outdir, & 282 prefix, wfcdir, max_seconds, ekin_conv_thr, etot_conv_thr, & 283 forc_conv_thr, pseudo_dir, disk_io, tefield, dipfield, lberry, & 284 gdir, nppstr, wf_collect, lelfield, nberrycyc, refg, & 285 tefield2, saverho, tabps, lkpoint_dir, use_wannier, lecrpa, & 286 tqmmm, vdw_table_name, lorbm, memory, point_label_type, & 287 lfcpopt, lfcpdyn, input_xml_schema_file, gate 288! 289!=----------------------------------------------------------------------------=! 290! SYSTEM Namelist Input Parameters 291!=----------------------------------------------------------------------------=! 292! 293 INTEGER :: ibrav = 14 294 ! index of the the Bravais lattice 295 ! Note: in variable cell CP molecular dynamics, usually one does 296 ! not want to put constraints on the cell symmetries, thus 297 ! ibrav = 14 is used 298 299 REAL(DP) :: celldm(6) = 0.0_DP 300 ! dimensions of the cell (lattice parameters and angles) 301 302 REAL(DP) :: a = 0.0_DP 303 REAL(DP) :: c = 0.0_DP 304 REAL(DP) :: b = 0.0_DP 305 REAL(DP) :: cosab = 0.0_DP 306 REAL(DP) :: cosac = 0.0_DP 307 REAL(DP) :: cosbc = 0.0_DP 308 ! Alternate definition of the cell - use either this or celldm 309 310 REAL(DP) :: ref_alat = 0.0_DP 311 ! reference cell alat in a.u. (see REF_CELL_PARAMETERS card) 312 313 INTEGER :: nat = 0 314 ! total number of atoms 315 316 INTEGER :: ntyp = 0 317 ! number of atomic species 318 319 INTEGER :: nbnd = 0 320 ! number of electronic states, this parameter is MANDATORY in CP 321 322 REAL(DP):: tot_charge = 0.0_DP 323 ! total system charge 324 325 REAL(DP) :: tot_magnetization = -1.0_DP 326 ! majority - minority spin. 327 ! A value < 0 means unspecified 328 329 REAL(DP) :: ecutwfc = 0.0_DP 330 ! energy cutoff for wave functions in k-space ( in Rydberg ) 331 ! this parameter is MANDATORY 332 333 REAL(DP) :: ecutrho = 0.0_DP 334 ! energy cutoff for charge density in k-space ( in Rydberg ) 335 ! by default its value is "4 * ecutwfc" 336 337 INTEGER :: nr1 = 0 338 INTEGER :: nr2 = 0 339 INTEGER :: nr3 = 0 340 ! dimensions of the real space grid for charge and potentials 341 ! presently NOT used in CP 342 343 INTEGER :: nr1s = 0 344 INTEGER :: nr2s = 0 345 INTEGER :: nr3s = 0 346 ! dimensions of the real space grid for wavefunctions 347 ! presently NOT used in CP 348 349 INTEGER :: nr1b = 0 350 INTEGER :: nr2b = 0 351 INTEGER :: nr3b = 0 352 ! dimensions of the "box" grid for Ultrasoft pseudopotentials 353 354 CHARACTER(len=80) :: occupations = 'fixed' 355 ! select the way electronic states are filled 356 ! See card 'OCCUPATIONS' if ocupations='from_input' 357 358 CHARACTER(len=80) :: smearing = 'gaussian' 359 ! select the way electronic states are filled for metalic systems 360 361 REAL(DP) :: degauss = 0.0_DP 362 ! parameter for the smearing functions - NOT used in CP 363 364 INTEGER :: nspin = 1 365 ! number of spinors 366 ! "nspin = 1" for LDA simulations 367 ! "nspin = 2" for LSD simulations 368 ! "nspin = 4" for NON COLLINEAR simulations 369 370 371 LOGICAL :: nosym = .true., noinv = .false. 372 ! (do not) use symmetry, q => -q symmetry in k-point generation 373 LOGICAL :: nosym_evc = .false. 374 ! if .true. use symmetry only to symmetrize k points 375 376 LOGICAL :: force_symmorphic = .false. 377 ! if .true. disable fractionary translations (nonsymmorphic groups) 378 LOGICAL :: use_all_frac = .false. 379 ! if .true. enable usage of all fractionary translations, 380 ! disabling check if they are commensurate with FFT grid 381 382 REAL(DP) :: ecfixed = 0.0_DP, qcutz = 0.0_DP, q2sigma = 0.0_DP 383 ! parameters for modified kinetic energy functional to be used 384 ! in variable-cell constant cut-off simulations 385 386 CHARACTER(len=80) :: input_dft = 'none' 387 ! Variable used to overwrite dft definition contained in 388 ! pseudopotential files; 'none' means DFT is read from pseudos. 389 ! Only used in PW - allowed values: any legal DFT value 390 391 REAL(DP) :: starting_charge( nsx ) = 0.0_DP 392 ! ONLY PW 393 394 REAL(DP) :: starting_magnetization( nsx ) = 0.0_DP 395 ! ONLY PW 396 397 LOGICAL :: lda_plus_u = .false. 398 ! Use DFT+U(+V) method - following are the needed parameters 399 INTEGER :: lda_plus_u_kind = 0 400 INTEGER :: lback(nsx) = -1 401 INTEGER :: l1back(nsx) = -1 402 INTEGER, PARAMETER :: nspinx=2, lqmax=7 403 REAL(DP) :: starting_ns_eigenvalue(lqmax,nspinx,nsx) = -1.0_DP 404 REAL(DP) :: hubbard_u(nsx) = 0.0_DP 405 REAL(DP) :: hubbard_u_back(nsx) = 0.0_DP 406 REAL(DP) :: hubbard_v(natx,natx*(2*sc_size+1)**3,4) = 0.0_DP 407 REAL(DP) :: hubbard_j0(nsx) = 0.0_DP 408 REAL(DP) :: hubbard_j(3,nsx) = 0.0_DP 409 REAL(DP) :: hubbard_alpha(nsx) = 0.0_DP 410 REAL(DP) :: hubbard_alpha_back(nsx) = 0.0_DP 411 REAL(DP) :: hubbard_beta(nsx) = 0.0_DP 412 CHARACTER(len=80) :: U_projection_type = 'atomic' 413 CHARACTER(len=80) :: Hubbard_parameters = 'input' 414 LOGICAL :: reserv(nsx) = .FALSE. 415 LOGICAL :: reserv_back(nsx) = .FALSE. 416 LOGICAL :: hub_pot_fix = .FALSE. 417 LOGICAL :: backall(nsx) = .FALSE. 418 419 LOGICAL :: la2F = .false. 420 ! For electron-phonon calculations 421 ! 422 LOGICAL :: step_pen=.false. 423 REAL(DP) :: A_pen(10,nspinx) = 0.0_DP 424 REAL(DP) :: sigma_pen(10) = 0.01_DP 425 REAL(DP) :: alpha_pen(10) = 0.0_DP 426 427 ! next group of variables PWSCF ONLY 428 ! 429 ! 430 REAL(DP) :: exx_fraction = -1.0_DP ! if negative, use defaults 431 REAL(DP) :: screening_parameter = -1.0_DP 432 INTEGER :: nqx1 = 0, nqx2 = 0, nqx3=0 ! use the same values as nk1, nk2, nk3 433 !gau-pbe in 434 REAL(DP) :: gau_parameter = -1.0_DP 435 !gau-pbe out 436 ! 437 CHARACTER(len=80) :: exxdiv_treatment = 'gygi-baldereschi' 438 ! define how ro cure the Coulomb divergence in EXX 439 ! Allowed values are: 440 CHARACTER(len=80) :: exxdiv_treatment_allowed(6) 441 DATA exxdiv_treatment_allowed / 'gygi-baldereschi', 'gygi-bald', 'g-b',& 442 'vcut_ws', 'vcut_spherical', 'none' / 443 ! 444 LOGICAL :: x_gamma_extrapolation = .true. 445 REAL(DP) :: yukawa = 0.0_DP 446 REAL(DP) :: ecutvcut = 0.0_DP 447 ! auxiliary variables to define exxdiv treatment 448 LOGICAL :: adaptive_thr = .FALSE. 449 REAL(DP) :: conv_thr_init = 0.001_DP 450 REAL(DP) :: conv_thr_multi = 0.1_DP 451 REAL(DP) :: ecutfock = -1.d0 452 ! variables used in Lin Lin's ACE and SCDM 453 REAL(DP) :: localization_thr = 0.0_dp, scdmden=1.0d0, scdmgrd=1.0d0 454 INTEGER :: nscdm = 1 455 INTEGER :: n_proj = 0 456 LOGICAL :: scdm=.FALSE. 457 LOGICAL :: ace=.TRUE. 458 459 ! parameters for external electric field 460 INTEGER :: edir = 0 461 REAL(DP) :: emaxpos = 0.0_DP 462 REAL(DP) :: eopreg = 0.0_DP 463 REAL(DP) :: eamp = 0.0_DP 464 465 ! TB parameters for charged plate representing the gate 466 ! and a possible potential barrier 467 ! added also below to the namelist 468 REAL(DP) :: zgate = 0.5 469 LOGICAL :: relaxz = .false. 470 LOGICAL :: block = .false. 471 REAL(DP) :: block_1 = 0.45 472 REAL(DP) :: block_2 = 0.55 473 REAL(DP) :: block_height = 0.1 474 475 ! Various parameters for noncollinear calculations 476 LOGICAL :: noncolin = .false. 477 LOGICAL :: lspinorb = .false. 478 LOGICAL :: lforcet=.FALSE. 479 LOGICAL :: starting_spin_angle=.FALSE. 480 REAL(DP) :: lambda = 1.0_DP 481 REAL(DP) :: fixed_magnetization(3) = 0.0_DP 482 REAL(DP) :: angle1(nsx) = 0.0_DP 483 REAL(DP) :: angle2(nsx) = 0.0_DP 484 INTEGER :: report =-1 485 LOGICAL :: no_t_rev = .FALSE. 486 487 CHARACTER(len=80) :: constrained_magnetization = 'none' 488 REAL(DP) :: B_field(3) = 0.0_DP 489 ! A fixed magnetic field defined by the vector B_field is added 490 ! to the exchange and correlation magnetic field. 491 492 CHARACTER(len=80) :: sic = 'none' 493 ! CP only - SIC correction (D'avezac Mauri) 494 ! Parameters for SIC calculation 495 REAL(DP) :: sic_epsilon = 0.0_DP 496 REAL(DP) :: sic_alpha = 0.0_DP 497 LOGICAL :: force_pairing = .false. 498 499 LOGICAL :: spline_ps = .false. 500 ! use spline interpolation for pseudopotential 501 LOGICAL :: one_atom_occupations=.false. 502 503 CHARACTER(len=80) :: assume_isolated = 'none' 504 ! possible corrections for isolated systems: 505 ! 'none', 'makov-payne', 'martyna-tuckerman', 'esm' 506 ! plus ENVIRON-specific: 507 ! 'slabx', 'slaby', 'slabz', 'pcc' 508 509 CHARACTER(len=80) :: vdw_corr = 'none' 510 ! semi-empirical van der Waals corrections 511 ! (not to be confused with nonlocal functionals, 512 ! specified in input_dft!). Default is 'none', allowed values: 513 ! 'dft-d' or 'grimme-d2' [S.Grimme, J.Comp.Chem. 27, 1787 (2006)] 514 ! 'ts', 'ts-vdW', 'tkatchenko-scheffler' 515 ! (Tkatchenko & Scheffler, Phys. Rev. Lett. 102, 073005 (2009)) 516 ! 'xdm' (Otero de la Roza and Johnson, J. Chem. Phys. 136 (2012) 174109) 517 518 LOGICAL :: london = .false. 519 ! OBSOLESCENT: same as vdw_corr='grimme-d2' 520 ! other DFT-D parameters ( see Modules/mm_dispersion.f90 ) 521 ! london_s6 = default global scaling parameter for PBE 522 ! london_c6 = user specified atomic C6 coefficients 523 ! london_rvdw = user specified atomic vdw radii 524 REAL ( DP ) :: london_s6 = 0.75_DP ,& 525 london_rcut = 200.00_DP , & 526 london_c6( nsx ) = -1.0_DP, & 527 london_rvdw( nsx ) = -1.0_DP 528 529 ! Grimme-D3 (DFT-D3) dispersion correction. 530 ! version=3 is Grimme-D3, version=2 is Grimme-D2, 531 ! version=3 dftd3_threebody = .false. similar to grimme d2 but with 532 ! the two_body parameters of d3 533 integer :: dftd3_version = 3 534 logical :: dftd3_threebody = .true. 535 536 LOGICAL :: ts_vdw = .false. 537 ! OBSOLESCENT: same as vdw_corr='Tkatchenko-Scheffler' 538 LOGICAL :: ts_vdw_isolated = .FALSE. 539 ! if .TRUE., TS-vdW correction for isolated system 540 ! if .FALSE., TS-vdW correction for periodic system 541 REAL(DP) :: ts_vdw_econv_thr = 1.0E-6_DP 542 ! convergence criterion for TS-vdW energy for periodic system 543 ! 544 545 LOGICAL :: xdm = .FALSE. 546 ! OBSOLESCENT: same as vdw_corr='xdm' 547 REAL(DP) :: xdm_a1 = 0.0_DP 548 REAL(DP) :: xdm_a2 = 0.0_DP 549 ! 550 CHARACTER(LEN=3) :: esm_bc = 'pbc' 551 ! 'pbc': ordinary calculation with periodic boundary conditions 552 ! 'bc1': vacuum-slab-vacuum 553 ! 'bc2': metal-slab-metal 554 ! 'bc3': vacuum-slab-metal 555 556 REAL(DP) :: esm_efield = 0.0_DP 557 ! applied electronic field [Ryd/a.u.] (used only for esm_bc='bc2') 558 559 REAL(DP) :: esm_w = 0.0_DP 560 ! position of effective screening medium from z0=L_z/2 [a.u.] 561 ! note: z1 is given by z1=z0+abs(esm_w) 562 563 REAL(DP) :: esm_a = 0.0_DP 564 ! smoothness parameter for smooth-ESM (exp(2a(z-z1))) 565 566 REAL(DP) :: esm_zb = 0.0_DP 567 ! smearing width for Ewald summation (RSUM) in smooth-ESM 568 569 INTEGER :: esm_nfit = 4 570 ! number of z-grid points for polynomial fitting at cell edge 571 572 LOGICAL :: esm_debug = .FALSE. 573 ! used to enable debug mode (output v_hartree and v_local) 574 575 INTEGER :: esm_debug_gpmax = 0 576 ! if esm_debug is .TRUE., calculate v_hartree and v_local 577 ! for abs(gp)<=esm_debug_gpmax (gp is integer and has tpiba unit) 578 579 REAL(DP) :: fcp_mu = 0.0_DP 580 ! target Fermi energy 581 582 REAL(DP) :: fcp_mass = 10000.0_DP 583 ! mass for the FCP 584 585 REAL(DP) :: fcp_tempw = 300.0_DP 586 ! target temperature for the FCP dynamics 587 588 CHARACTER(LEN=8) :: fcp_relax = 'lm' 589 ! 'lm': line minimisation 590 ! 'mdiis': MDIIS algorithm 591 592 CHARACTER(len=8) :: fcp_relax_allowed(2) 593 DATA fcp_relax_allowed / 'lm', 'mdiis' / 594 595 REAL(DP) :: fcp_relax_step = 0.5_DP 596 ! step size for steepest descent 597 598 REAL(DP) :: fcp_relax_crit = 0.001_DP 599 ! threshold for force acting on FCP 600 601 INTEGER :: fcp_mdiis_size = 4 602 ! size of MDIIS algorism 603 604 REAL(DP) :: fcp_mdiis_step = 0.2_DP 605 ! step width of MDIIS algorism 606 607 INTEGER :: space_group = 0 608 ! space group number for coordinates given in crystallographic form 609 ! 610 LOGICAL :: uniqueb=.FALSE. 611 ! if .TRUE. for monoclinic lattice choose the b unique primitive 612 ! vectors 613 ! 614 INTEGER :: origin_choice = 1 615 ! for space groups that have more than one origin choice, choose 616 ! the origin (can be 1 or 2) 617 ! 618 LOGICAL :: rhombohedral = .TRUE. 619 ! 620 ! if .TRUE. for rhombohedral space groups give the coordinates 621 ! in rhombohedral axes. If .FALSE. in hexagonal axes, that are 622 ! converted internally in rhombohedral axes. 623 ! 624 625 626 NAMELIST / system / ibrav, celldm, a, b, c, cosab, cosac, cosbc, nat, & 627 ntyp, nbnd, ecutwfc, ecutrho, nr1, nr2, nr3, nr1s, nr2s, & 628 nr3s, nr1b, nr2b, nr3b, nosym, nosym_evc, noinv, use_all_frac, & 629 force_symmorphic, starting_charge, starting_magnetization, & 630 occupations, degauss, nspin, ecfixed, & 631 qcutz, q2sigma, lda_plus_U, lda_plus_u_kind, & 632 Hubbard_U, Hubbard_U_back, Hubbard_J, Hubbard_alpha, & 633 Hubbard_alpha_back, Hubbard_J0, Hubbard_beta, & 634 hub_pot_fix, Hubbard_V, Hubbard_parameters, & 635 backall, lback, l1back, reserv, reserv_back, & 636 edir, emaxpos, eopreg, eamp, smearing, starting_ns_eigenvalue, & 637 U_projection_type, input_dft, la2F, assume_isolated, & 638 nqx1, nqx2, nqx3, ecutfock, localization_thr, scdm, ace, & 639 scdmden, scdmgrd, nscdm, n_proj, & 640 exxdiv_treatment, x_gamma_extrapolation, yukawa, ecutvcut, & 641 exx_fraction, screening_parameter, ref_alat, & 642 noncolin, lspinorb, starting_spin_angle, lambda, angle1, angle2, & 643 report, lforcet, & 644 constrained_magnetization, B_field, fixed_magnetization, & 645 sic, sic_epsilon, force_pairing, sic_alpha, & 646 tot_charge, tot_magnetization, spline_ps, one_atom_occupations, & 647 vdw_corr, london, london_s6, london_rcut, london_c6, london_rvdw,& 648 dftd3_version, dftd3_threebody, & 649 ts_vdw, ts_vdw_isolated, ts_vdw_econv_thr, & 650 xdm, xdm_a1, xdm_a2, & 651 step_pen, A_pen, sigma_pen, alpha_pen, no_t_rev, & 652 esm_bc, esm_efield, esm_w, esm_nfit, esm_debug, esm_debug_gpmax, & 653 esm_a, esm_zb, fcp_mu, fcp_mass, fcp_tempw, fcp_relax, & 654 fcp_relax_step, fcp_relax_crit, fcp_mdiis_size, fcp_mdiis_step, & 655 space_group, uniqueb, origin_choice, rhombohedral, & 656 zgate, relaxz, block, block_1, block_2, block_height 657 658!=----------------------------------------------------------------------------=! 659! ELECTRONS Namelist Input Parameters 660!=----------------------------------------------------------------------------=! 661 662 REAL(DP) :: emass = 0.0_DP 663 ! effective electron mass in the CP Lagrangian, 664 ! in atomic units ( 1 a.u. of mass = 1/1822.9 a.m.u. = 9.10939 * 10^-31 kg ) 665 ! Typical values in CP simulation are between 100. and 1000. 666 667 REAL(DP) :: emass_cutoff = 0.0_DP 668 ! mass cut-off (in Rydbergs) for the Fourier acceleration 669 ! effective mass is rescaled for "G" vector components with kinetic 670 ! energy above "emass_cutoff" 671 ! Use a value grether than "ecutwfc" to disable Fourier acceleration. 672 673 CHARACTER(len=80) :: orthogonalization = 'ortho' 674 ! orthogonalization = 'Gram-Schmidt' | 'ortho'* 675 ! selects the orthonormalization method for electronic wave functions 676 ! 'Gram-Schmidt' use Gram-Schmidt algorithm 677 ! 'ortho' use iterative algorithm 678 679 REAL(DP) :: ortho_eps = 1.E-9_DP 680 ! meaningful only if orthogonalization = 'ortho' 681 ! tolerance for iterative orthonormalization, 682 ! a value of 1.d-8 is usually sufficent 683 684 INTEGER :: ortho_max = 50 685 ! meaningful only if orthogonalization = 'ortho' 686 ! maximum number of iterations for orthonormalization 687 ! usually between 20 and 300. 688 689 INTEGER :: electron_maxstep = 1000 690 ! maximum number of steps in electronic minimization 691 ! This parameter apply only when using 'cg' electronic or 692 ! ionic dynamics and electron_dynamics = 'CP-BO' 693 LOGICAL :: scf_must_converge = .true. 694 ! stop or continue if SCF does not converge 695 696 CHARACTER(len=80) :: electron_dynamics = 'none' 697 ! set how electrons should be moved 698 CHARACTER(len=80) :: electron_dynamics_allowed(7) 699 DATA electron_dynamics_allowed & 700 / 'default', 'sd', 'cg', 'damp', 'verlet', 'none', 'cp-bo' / 701 702 REAL(DP) :: electron_damping = 0.0_DP 703 ! meaningful only if " electron_dynamics = 'damp' " 704 ! damping frequency times delta t, optimal values could be 705 ! calculated with the formula 706 ! sqrt(0.5*log((E1-E2)/(E2-E3))) 707 ! where E1 E2 E3 are successive values of the DFT total energy 708 ! in a steepest descent simulations 709 710 CHARACTER(len=80) :: electron_velocities = 'default' 711 ! electron_velocities = 'zero' | 'default'* 712 ! 'zero' restart setting electronic velocities to zero 713 ! 'default' restart using electronic velocities of the previous run 714 715 CHARACTER(len=80) :: electron_temperature = 'not_controlled' 716 ! electron_temperature = 'nose' | 'not_controlled'* | 'rescaling' 717 ! 'nose' control electronic temperature using Nose thermostat 718 ! see parameter "fnosee" and "ekincw" 719 ! 'rescaling' control electronic temperature via velocities rescaling 720 ! 'not_controlled' electronic temperature is not controlled 721 722 REAL(DP) :: ekincw = 0.0_DP 723 ! meaningful only with "electron_temperature /= 'not_controlled' " 724 ! value of the average kinetic energy (in atomic units) forced 725 ! by the temperature control 726 727 REAL(DP) :: fnosee = 0.0_DP 728 ! meaningful only with "electron_temperature = 'nose' " 729 ! oscillation frequency of the nose thermostat (in terahertz) 730 731 CHARACTER(len=80) :: startingwfc = 'random' 732 ! startingwfc = 'atomic' | 'atomic+random' | 'random' | 'file' 733 ! define how the code should initialize the wave function 734 ! 'atomic' start from superposition of atomic wave functions 735 ! 'atomic+random' as above, plus randomization 736 ! 'random' start from random wave functions 737 ! 'file' read wavefunctions from file 738 739 REAL(DP) :: ampre = 0.0_DP 740 ! meaningful only if "startingwfc = 'random'", amplitude of the 741 ! randomization ( allowed values: 0.0 - 1.0 ) 742 743 REAL(DP) :: grease = 0.0_DP 744 ! a number <= 1, very close to 1: the damping in electronic 745 ! damped dynamics is multiplied at each time step by "grease" 746 ! (avoids overdamping close to convergence: Obsolete ?) 747 ! grease = 1 : normal damped dynamics 748 ! used only in CP 749 750 INTEGER :: diis_size = 0 751 ! meaningful only with " electron_dynamics = 'diis' " 752 ! size of the matrix used for the inversion in the iterative subspace 753 ! default is 4, allowed value 1-5 754 755 INTEGER :: diis_nreset = 0 756 ! meaningful only with " electron_dynamics = 'diis' " 757 ! number of steepest descendent step after a reset of the diis 758 ! iteration, default value is 3 759 760 REAL(DP) :: diis_hcut = 0.0_DP 761 ! meaningful only with " electron_dynamics = 'diis' " 762 ! energy cutoff (a.u.), above which an approximate diagonal 763 ! hamiltonian is used in finding the direction to the minimum 764 ! default is "1.0" 765 766 REAL(DP) :: diis_wthr = 1.E-4_DP 767 ! meaningful only with " electron_dynamics = 'diis' " 768 ! convergence threshold for wave function 769 ! this criterion is satisfied when the maximum change 770 ! in the wave functions component between two diis steps 771 ! is less than this threshold 772 ! default value is ekin_conv_thr 773 774 REAL(DP) :: diis_delt = 1.0_DP 775 ! meaningful only with " electron_dynamics = 'diis' " 776 ! electronic time step used in the steepest descendent step 777 ! default is "dt" 778 779 INTEGER :: diis_maxstep = 100 780 ! meaningful only with " electron_dynamics = 'diis' " 781 ! maximum number of iteration in the diis minimization 782 ! default is electron_maxstep 783 784 LOGICAL :: diis_rot = .false. 785 ! meaningful only with " electron_dynamics = 'diis' " 786 ! if "diis_rot = .TRUE." enable diis with charge mixing and rotations 787 ! default is "diis_rot = .FALSE." 788 789 REAL(DP) :: diis_fthr = 1.E-3_DP 790 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 791 ! convergence threshold for ionic force 792 ! this criterion is satisfied when the maximum change 793 ! in the atomic force between two diis steps 794 ! is less than this threshold 795 ! default value is "0.0" 796 797 REAL(DP) :: diis_temp = 0.0_DP 798 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 799 ! electronic temperature, significant only if ??? 800 801 REAL(DP) :: diis_achmix = 0.0_DP 802 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 803 ! "A" parameter in the charge mixing formula 804 ! chmix = A * G^2 / (G^2 + G0^2) , G represents reciprocal lattice vectors 805 806 REAL(DP) :: diis_g0chmix = 0.0_DP 807 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 808 ! "G0^2" parameter in the charge mixing formula 809 810 INTEGER :: diis_nchmix = 0 811 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 812 ! dimension of the charge mixing 813 814 REAL(DP) :: diis_g1chmix = 0.0_DP 815 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 816 ! "G1^2" parameter in the charge mixing formula 817 ! metric = (G^2 + G1^2) / G^2 , G represents reciprocal lattice vectors 818 819 INTEGER :: diis_nrot(3) = 0 820 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 821 ! start upgrading the charge density every "diis_nrot(1)" steps, 822 ! then every "diis_nrot(2)", and at the end every "diis_nrot(3)", 823 ! depending on "diis_rothr" 824 825 REAL(DP) :: diis_rothr(3) = 1.E-4_DP 826 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 827 ! threshold on the charge difference between two diis step 828 ! when max charge difference is less than "diis_rothr(1)", switch 829 ! between the "diis_nrot(1)" upgrade frequency to "diis_nrot(2)", 830 ! then when the max charge difference is less than "diis_rothr(2)", 831 ! switch between "diis_nrot(2)" and "diis_nrot(3)", upgrade frequency, 832 ! finally when the max charge difference is less than "diis_nrot(3)" 833 ! convergence is achieved 834 835 REAL(DP) :: diis_ethr = 1.E-4_DP 836 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 837 ! convergence threshold for energy 838 ! this criterion is satisfied when the change 839 ! in the energy between two diis steps 840 ! is less than this threshold 841 ! default value is etot_conv_thr 842 843 LOGICAL :: diis_chguess = .false. 844 ! meaningful only with "electron_dynamics='diis' " and "diis_rot=.TRUE." 845 ! if "diis_chguess = .TRUE." enable charge density guess 846 ! between two diis step, defaut value is "diis_chguess = .FALSE." 847 848 CHARACTER(len=80) :: mixing_mode = 'default' 849 ! type of mixing algorithm for charge self-consistency 850 ! used only in PWscf 851 852 REAL(DP) :: mixing_beta = 0.0_DP 853 ! parameter for mixing algorithm 854 ! used only in PWscf 855 856 INTEGER :: mixing_ndim = 0 857 ! dimension of mixing subspace 858 ! used only in PWscf 859 860 CHARACTER(len=80) :: diagonalization = 'david' 861 ! diagonalization = 'david', 'cg' or 'ppcg' 862 ! algorithm used by PWscf for iterative diagonalization 863 864 REAL(DP) :: diago_thr_init = 0.0_DP 865 ! convergence threshold for the first iterative diagonalization. 866 ! used only in PWscf 867 868 INTEGER :: diago_cg_maxiter = 100 869 ! max number of iterations for the first iterative diagonalization 870 ! using conjugate-gradient algorithm - used only in PWscf 871 872 INTEGER :: diago_ppcg_maxiter = 100 873 ! max number of iterations for the first iterative diagonalization 874 ! using projected preconditioned conjugate-gradient algorithm - used only in PWscf 875 876 INTEGER :: diago_david_ndim = 4 877 ! dimension of the subspace used in Davidson diagonalization 878 ! used only in PWscf 879 880 LOGICAL :: diago_full_acc = .false. 881 882 REAL(DP) :: conv_thr = 1.E-6_DP 883 ! convergence threshold in electronic ONLY minimizations 884 ! used only in PWscf 885 ! used in electron_dynamics = 'CP-BO' in CP and PWscf 886 887 INTEGER :: mixing_fixed_ns = 0 888 ! For DFT+U calculations, PWscf only 889 890 CHARACTER(len=80) :: startingpot = 'potfile' 891 ! specify the file containing the DFT potential of the system 892 ! used only in PWscf 893 894 INTEGER :: n_inner = 2 895 ! number of inner loop per CG iteration. 896 ! used only in CP 897 898 INTEGER :: niter_cold_restart = 1 899 !frequency of full cold smearing inner cycle (in iterations) 900 901 REAL(DP) :: lambda_cold 902 !step for not complete cold smearing inner cycle 903 904 LOGICAL :: tgrand = .false. 905 ! whether to do grand-canonical calculations. 906 907 REAL(DP) :: fermi_energy = 0.0_DP 908 ! chemical potential of the grand-canonical ensemble. 909 910 CHARACTER(len=80) :: rotation_dynamics = "line-minimization" 911 ! evolution the rotational degrees of freedom. 912 913 CHARACTER(len=80) :: occupation_dynamics = "line-minimization" 914 ! evolution of the occupational degrees of freedom. 915 916 REAL(DP) :: rotmass = 0 917 ! mass for the rotational degrees of freedom. 918 919 REAL(DP) :: occmass = 0 920 ! mass for the occupational degrees of freedom. 921 922 REAL(DP) :: occupation_damping = 0 923 ! damping for the rotational degrees of freedom. 924 925 REAL(DP) :: rotation_damping = 0 926 ! damping for the occupational degrees of freedom. 927 928 LOGICAL :: tcg = .true. 929 ! if true perform in cpv conjugate gradient minimization of electron energy 930 931 INTEGER :: maxiter = 100 932 ! max number of conjugate gradient iterations 933 934 REAL(DP) :: etresh =1.0E-7_DP 935 ! treshhold on energy 936 937 REAL(DP) :: passop =0.3_DP 938 ! small step for parabolic interpolation 939 940 INTEGER :: niter_cg_restart 941 !frequency of restart for the conjugate gradient algorithm in iterations 942 943 INTEGER :: epol = 3 944 ! electric field direction 945 946 REAL(DP) :: efield =0.0_DP 947 ! electric field intensity in atomic units 948 949 ! real_space routines for US pps 950 LOGICAL :: real_space = .false. 951 952 953 REAL(DP) :: efield_cart(3) 954 ! electric field vector in cartesian system of reference 955 956 CHARACTER(len=80) :: efield_phase='none' 957 ! for Berry's phase electric field selection of string phases 958 959 INTEGER :: epol2 = 3 960 ! electric field direction 961 962 REAL(DP) :: efield2 =0.0_DP 963 ! electric field intensity in atomic units 964 965 LOGICAL :: tqr = .false. 966 ! US contributions are added in real space 967 968 LOGICAL :: tq_smoothing = .false. 969 ! US augmentation charge is smoothed before use 970 971 LOGICAL :: tbeta_smoothing = .false. 972 ! beta function are smoothed before use 973 974 LOGICAL :: occupation_constraints = .false. 975 ! If true perform CP dynamics with constrained occupations 976 ! to be used together with penalty functional ... 977 978 ! 979 ! ... CP-BO ... 980 LOGICAL :: tcpbo = .FALSE. 981 ! if true perform CP-BO minimization of electron energy 982 983 REAL(DP) :: emass_emin = 0.0_DP 984 ! meaningful only if electron_dynamics = 'CP-BO' 985 ! effective electron mass used in CP-BO electron minimization in 986 ! atomic units ( 1 a.u. of mass = 1/1822.9 a.m.u. = 9.10939 * 10^-31 kg ) 987 988 REAL(DP) :: emass_cutoff_emin = 0.0_DP 989 ! meaningful only if electron_dynamics = 'CP-BO' 990 ! mass cut-off (in Rydbergs) for the Fourier acceleration in CP-BO 991 ! electron minimization 992 993 REAL(DP) :: electron_damping_emin = 0.0_DP 994 ! meaningful only if electron_dynamics = 'CP-BO' 995 ! damping parameter utilized in CP-BO electron minimization 996 997 REAL(DP) :: dt_emin = 0.0_DP 998 ! meaningful only if electron_dynamics = 'CP-BO' 999 ! time step for CP-BO electron minimization dynamics, in atomic units 1000 ! CP: 1 a.u. of time = 2.4189 * 10^-17 s, PW: twice that much 1001 1002 NAMELIST / electrons / emass, emass_cutoff, orthogonalization, & 1003 electron_maxstep, scf_must_converge, ortho_eps, ortho_max, electron_dynamics, & 1004 electron_damping, electron_velocities, electron_temperature, & 1005 ekincw, fnosee, ampre, grease, & 1006 diis_size, diis_nreset, diis_hcut, & 1007 diis_wthr, diis_delt, diis_maxstep, diis_rot, diis_fthr, & 1008 diis_temp, diis_achmix, diis_g0chmix, diis_g1chmix, & 1009 diis_nchmix, diis_nrot, diis_rothr, diis_ethr, diis_chguess, & 1010 mixing_mode, mixing_beta, mixing_ndim, mixing_fixed_ns, & 1011 tqr, tq_smoothing, tbeta_smoothing, & 1012 diago_cg_maxiter, diago_david_ndim, diagonalization, & 1013 startingpot, startingwfc , conv_thr, & 1014 adaptive_thr, conv_thr_init, conv_thr_multi, & 1015 diago_thr_init, n_inner, fermi_energy, rotmass, occmass, & 1016 rotation_damping, occupation_damping, rotation_dynamics, & 1017 occupation_dynamics, tcg, maxiter, etresh, passop, epol, & 1018 efield, epol2, efield2, diago_full_acc, & 1019 occupation_constraints, niter_cg_restart, & 1020 niter_cold_restart, lambda_cold, efield_cart, real_space, & 1021 tcpbo,emass_emin, emass_cutoff_emin, electron_damping_emin, & 1022 dt_emin, efield_phase 1023 1024! 1025!=----------------------------------------------------------------------------=! 1026! IONS Namelist Input Parameters 1027!=----------------------------------------------------------------------------=! 1028 1029 1030 CHARACTER(len=80) :: ion_dynamics = 'none' 1031 ! set how ions should be moved 1032 CHARACTER(len=80) :: ion_dynamics_allowed(10) 1033 DATA ion_dynamics_allowed / 'none', 'sd', 'cg', 'langevin', & 1034 'damp', 'verlet', 'bfgs', 'beeman',& 1035 'langevin-smc', 'ipi' / 1036 1037 REAL(DP) :: ion_radius(nsx) = 0.5_DP 1038 ! pseudo-atomic radius of the i-th atomic species (CP only) 1039 ! for Ewald summation: typical values range between 0.5 and 2.0 1040 INTEGER :: iesr = 1 1041 ! perform Ewald summation on iesr*iesr*iesr cells - CP only 1042 1043 REAL(DP) :: ion_damping = 0.2_DP 1044 ! meaningful only if " ion_dynamics = 'damp' " 1045 ! damping frequency times delta t, optimal values could be 1046 ! calculated with the formula 1047 ! sqrt(0.5*log((E1-E2)/(E2-E3))) 1048 ! where E1 E2 E3 are successive values of the DFT total energy 1049 ! in a ionic steepest descent simulation 1050 1051 CHARACTER(len=80) :: ion_positions = 'default' 1052 ! ion_positions = 'default'* | 'from_input' 1053 ! 'default' restart the simulation with atomic positions read 1054 ! from the restart file 1055 ! 'from_input' restart the simulation with atomic positions read 1056 ! from standard input ( see the card 'ATOMIC_POSITIONS' ) 1057 1058 CHARACTER(len=80) :: ion_velocities = 'default' 1059 ! ion_velocities = 'zero' | 'default'* | 'random' | 'from_input' 1060 ! 'default' restart the simulation with atomic velocities read 1061 ! from the restart file 1062 ! 'random' start the simulation with random atomic velocities 1063 ! 'from_input' restart the simulation with atomic velocities read 1064 ! from standard input (see the card 'ATOMIC_VELOCITIES' ) 1065 ! 'zero' restart the simulation with atomic velocities set to zero 1066 1067 CHARACTER(len=80) :: ion_temperature = 'not_controlled' 1068 ! ion_temperature = 'nose' | 'not_controlled'* | 'rescaling' | 1069 ! 'berendsen' | 'andersen' | 'rescale-v' | 'rescale-T' | 'reduce-T' 1070 ! 1071 ! 'nose' control ionic temperature using Nose thermostat 1072 ! see parameters "fnosep" and "tempw" 1073 ! 'rescaling' control ionic temperature via velocity rescaling 1074 ! see parameters "tempw" and "tolp" 1075 ! 'rescale-v' control ionic temperature via velocity rescaling 1076 ! see parameters "tempw" and "nraise" 1077 ! 'rescale-T' control ionic temperature via velocity rescaling 1078 ! see parameter "delta_t" 1079 ! 'reduce-T' reduce ionic temperature 1080 ! see parameters "nraise", delta_t" 1081 ! 'berendsen' control ionic temperature using "soft" velocity 1082 ! rescaling - see parameters "tempw" and "nraise" 1083 ! 'andersen' control ionic temperature using Andersen thermostat 1084 ! see parameters "tempw" and "nraise" 1085 ! 'not_controlled' ionic temperature is not controlled 1086 1087 REAL(DP) :: tempw = 300.0_DP 1088 ! meaningful only with "ion_temperature /= 'not_controlled' " 1089 ! value of the ionic temperature (in Kelvin) forced 1090 ! by the temperature control 1091 1092 INTEGER, PARAMETER :: nhclm = 4 1093 REAL(DP) :: fnosep( nhclm ) = 50.0_DP 1094 ! meaningful only with "ion_temperature = 'nose' " 1095 ! oscillation frequency of the nose thermostat (in terahertz) 1096 ! nhclm is the max length for the chain; it can be easily increased 1097 ! since the restart file should be able to handle it 1098 ! perhaps better to align nhclm by 4 1099 1100 INTEGER :: nhpcl = 0 1101 ! non-zero only with "ion_temperature = 'nose' " 1102 ! this defines the length of the Nose-Hoover chain 1103 1104 INTEGER :: nhptyp = 0 1105 ! this parameter set the nose hoover thermostat to more than one 1106 1107 INTEGER :: nhgrp(nsx)=0 1108 ! this is the array to assign thermostats to atomic types 1109 ! allows to use various thermostat setups 1110 1111 INTEGER :: ndega = 0 1112 ! this is the parameter to control active degrees of freedom 1113 ! used for temperature control and the Nose-Hoover chains 1114 1115 REAL(DP) :: tolp = 50.0_DP 1116 ! meaningful only with "ion_temperature = 'rescaling' " 1117 ! tolerance (in Kelvin) of the rescaling. When ionic temperature 1118 ! differs from "tempw" more than "tolp" apply rescaling. 1119 1120 REAL(DP) :: fnhscl(nsx)=-1.0_DP 1121 ! this is to scale the target energy, in case there are constraints 1122 ! the dimension is the same as nhgrp, meaning that atomic type 1123 ! i with a group nhgrp(i) is scaled by fnhscl(i) 1124 1125 LOGICAL :: tranp(nsx) = .false. 1126 ! tranp(i) control the randomization of the i-th atomic specie 1127 ! .TRUE. randomize ionic positions ( see "amprp" ) 1128 ! .FALSE. do nothing 1129 1130 REAL(DP) :: amprp(nsx) = 0.0_DP 1131 ! amprp(i) meaningful only if "tranp(i) = .TRUE.", amplitude of the 1132 ! randomization ( allowed values: 0.0 - 1.0 ) for the i-th atomic specie. 1133 ! Add to the positions a random displacements vector ( in bohr radius ) 1134 ! defined as: amprp( i ) * ( X, Y, Z ) 1135 ! where X, Y, Z are pseudo random number in the interval [ -0.5 , 0.5 ] 1136 1137 REAL(DP) :: greasp = 0.0_DP 1138 ! same as "grease", for ionic damped dynamics 1139 ! NOT used in FPMD 1140 1141 INTEGER :: ion_nstepe = 1 1142 ! number of electronic steps for each ionic step 1143 1144 INTEGER :: ion_maxstep = 1000 1145 ! maximum number of step in ionic minimization 1146 1147 REAL(DP) :: upscale = 100.0_DP 1148 ! Max reduction allowed in scf threshold during optimization 1149 1150 CHARACTER(len=80) :: pot_extrapolation = 'default', & 1151 wfc_extrapolation = 'default' 1152 ! These variables are used only by PWSCF 1153 1154 LOGICAL :: refold_pos 1155 LOGICAL :: remove_rigid_rot = .false. 1156 1157 ! 1158 ! ... delta_T, nraise, tolp are used to change temperature in PWscf 1159 ! 1160 1161 REAL(DP) :: delta_t = 1.0_DP 1162 1163 INTEGER :: nraise = 1 1164 1165 ! 1166 ! ... variables added for new BFGS algorithm 1167 ! 1168 1169 INTEGER :: bfgs_ndim = 1 1170 1171 REAL(DP) :: trust_radius_max = 0.8_DP 1172 REAL(DP) :: trust_radius_min = 1.E-3_DP 1173 REAL(DP) :: trust_radius_ini = 0.5_DP 1174 1175 REAL(DP) :: w_1 = 0.5E-1_DP 1176 REAL(DP) :: w_2 = 0.5_DP 1177 1178 LOGICAL :: l_mplathe=.false. !if true apply Muller Plathe strategy 1179 INTEGER :: n_muller=0!number of intermediate sub-cells 1180 INTEGER :: np_muller=1!period for velocity exchange 1181 LOGICAL :: l_exit_muller=.false.!if true do muller exchange after last MD step 1182 1183 1184 ! 1185 NAMELIST / ions / ion_dynamics, iesr, ion_radius, ion_damping, & 1186 ion_positions, ion_velocities, ion_temperature, & 1187 tempw, fnosep, nhgrp, fnhscl, nhpcl, nhptyp, ndega, tranp, & 1188 amprp, greasp, tolp, ion_nstepe, ion_maxstep, & 1189 refold_pos, upscale, delta_t, pot_extrapolation, & 1190 wfc_extrapolation, nraise, remove_rigid_rot, & 1191 trust_radius_max, trust_radius_min, & 1192 trust_radius_ini, w_1, w_2, bfgs_ndim,l_mplathe, & 1193 n_muller,np_muller,l_exit_muller 1194 1195 1196!=----------------------------------------------------------------------------=! 1197! CELL Namelist Input Parameters 1198!=----------------------------------------------------------------------------=! 1199! 1200 CHARACTER(len=80) :: cell_parameters = 'default' 1201 ! cell_parameters = 'default'* | 'from_input' 1202 ! 'default' restart the simulation with cell parameters read 1203 ! from the restart file or "celldm" if 1204 ! "restart = 'from_scratch'" 1205 ! 'from_input' restart the simulation with cell parameters 1206 ! from standard input ( see the card 'CELL_PARAMETERS' ) 1207 1208 CHARACTER(len=80) :: cell_dynamics = 'none' 1209 ! set how the cell should be moved 1210 CHARACTER(len=80) :: cell_dynamics_allowed(8) 1211 DATA cell_dynamics_allowed / 'sd', 'pr', 'none', 'w', 'damp-pr', & 1212 'damp-w', 'bfgs', 'ipi' / 1213 1214 CHARACTER(len=80) :: cell_velocities = 'default' 1215 ! cell_velocities = 'zero' | 'default'* 1216 ! 'zero' restart setting cell velocitiy to zero 1217 ! 'default' restart using cell velocity of the previous run 1218 1219 REAL(DP) :: press = 0.0_DP 1220 ! external pressure (in GPa, remember 1 kbar = 10^8 Pa) 1221 1222 REAL(DP) :: wmass = 0.0_DP 1223 ! effective cell mass in the Parrinello-Rahman Lagrangian (in atomic units) 1224 ! of the order of magnitude of the total atomic mass 1225 ! (sum of the mass of the atoms) within the simulation cell. 1226 ! if you do not specify this parameters, the code will compute 1227 ! its value based on some physical consideration 1228 1229 CHARACTER(len=80) :: cell_temperature = 'not_controlled' 1230 ! cell_temperature = 'nose' | 'not_controlled'* | 'rescaling' 1231 ! 'nose' control cell temperature using Nose thermostat 1232 ! see parameters "fnoseh" and "temph" 1233 ! 'rescaling' control cell temperature via velocities rescaling 1234 ! 'not_controlled' cell temperature is not controlled 1235 ! NOT used in FPMD 1236 1237 REAL(DP) :: temph = 0.0_DP 1238 ! meaningful only with "cell_temperature /= 'not_controlled' " 1239 ! value of the cell temperature (in Kelvin) forced 1240 ! by the temperature control 1241 1242 REAL(DP) :: fnoseh = 1.0_DP 1243 ! meaningful only with "cell_temperature = 'nose' " 1244 ! oscillation frequency of the nose thermostat (in terahertz) 1245 1246 REAL(DP) :: greash = 0.0_DP 1247 ! same as "grease", for cell damped dynamics 1248 1249 CHARACTER(len=80) :: cell_dofree = 'all' 1250 ! cell_dofree = 'all'* | 'volume' | 'x' | 'y' | 'z' | 'xy' | 'xz' | 'yz' | 'xyz' 1251 ! select which of the cell parameters should be moved 1252 ! 'all' all axis and angles are propagated (default) 1253 ! 'volume' the cell is simply rescaled, without changing the shape 1254 ! 'x' only the "x" axis is moved 1255 ! 'y' only the "y" axis is moved 1256 ! 'z' only the "z" axis is moved 1257 ! 'xy' only the "x" and "y" axis are moved, angles are unchanged 1258 ! 'xz' only the "x" and "z" axis are moved, angles are unchanged 1259 ! 'yz' only the "y" and "z" axis are moved, angles are unchanged 1260 ! 'xyz' "x", "y" and "z" axis are moved, angles are unchanged 1261 1262 REAL(DP) :: cell_factor = 0.0_DP 1263 ! NOT used in FPMD 1264 1265 INTEGER :: cell_nstepe = 1 1266 ! number of electronic steps for each cell step 1267 1268 REAL(DP) :: cell_damping = 0.1_DP 1269 ! meaningful only if " cell_dynamics = 'damp' " 1270 ! damping frequency times delta t, optimal values could be 1271 ! calculated with the formula 1272 ! sqrt(0.5*log((E1-E2)/(E2-E3))) 1273 ! where E1 E2 E3 are successive values of the DFT total energy 1274 ! in a ionic steepest descent simulation 1275 1276 REAL(DP) :: press_conv_thr = 0.5_DP 1277 1278 LOGICAL :: treinit_gvecs = .FALSE. 1279 ! if true all the quantities related to fft g vectors are updated at 1280 ! step of variable cell structural optimization 1281 1282 NAMELIST / cell / cell_parameters, cell_dynamics, cell_velocities, & 1283 press, wmass, cell_temperature, temph, fnoseh, & 1284 cell_dofree, greash, cell_factor, cell_nstepe, & 1285 cell_damping, press_conv_thr, treinit_gvecs 1286 1287! 1288!=----------------------------------------------------------------------------=!! 1289! PRESS_AI Namelist Input Parameters 1290!=----------------------------------------------------------------------------=! 1291! 1292! 1293 LOGICAL :: abivol = .false. 1294 LOGICAL :: abisur = .false. 1295 LOGICAL :: pvar = .false. 1296 LOGICAL :: fill_vac=.false. 1297 LOGICAL :: scale_at=.false. 1298 LOGICAL :: t_gauss =.false. 1299 LOGICAL :: jellium= .false. 1300 LOGICAL :: cntr(nsx)=.false. 1301 REAL(DP) :: P_ext = 0.0_DP 1302 REAL(DP) :: P_in = 0.0_DP 1303 REAL(DP) :: P_fin = 0.0_DP 1304 REAL(DP) :: rho_thr = 0.0_DP 1305 REAL(DP) :: step_rad(nsx) = 0.0_DP 1306 REAL(DP) :: Surf_t = 0.0_DP 1307 REAL(DP) :: dthr = 0.0_DP 1308 REAL(DP) :: R_j = 0.0_DP 1309 REAL(DP) :: h_j = 0.0_DP 1310 REAL(DP) :: delta_eps = 0.0_DP 1311 REAL(DP) :: delta_sigma=0.0_DP 1312 INTEGER :: n_cntr = 0 1313 INTEGER :: axis = 0 1314 1315 NAMELIST / press_ai / abivol, P_ext, pvar, P_in, P_fin, rho_thr, & 1316 & step_rad, delta_eps, delta_sigma, n_cntr, & 1317 & fill_vac, scale_at, t_gauss, abisur, & 1318 & Surf_t, dthr, cntr, axis, jellium, R_j, h_j 1319 1320!=----------------------------------------------------------------------------=! 1321! WANNIER Namelist Input Parameters 1322!=----------------------------------------------------------------------------=! 1323 1324 LOGICAL :: wf_efield 1325 LOGICAL :: wf_switch 1326 ! 1327 INTEGER :: sw_len 1328 ! 1329 REAL(DP) :: efx0, efy0, efz0 1330 REAL(DP) :: efx1, efy1, efz1 1331 ! 1332 INTEGER :: wfsd 1333 ! 1334 REAL(DP) :: wfdt 1335 REAL(DP) :: maxwfdt 1336 REAL(DP) :: wf_q 1337 REAL(DP) :: wf_friction 1338!======================================================================= 1339!exx_wf related 1340 INTEGER :: vnbsp 1341 INTEGER :: exx_neigh 1342 REAL(DP) :: exx_poisson_eps 1343 REAL(DP) :: exx_dis_cutoff 1344 REAL(DP) :: exx_ps_rcut_self 1345 REAL(DP) :: exx_ps_rcut_pair 1346 REAL(DP) :: exx_me_rcut_self 1347 REAL(DP) :: exx_me_rcut_pair 1348 LOGICAL :: exx_use_cube_domain 1349!======================================================================= 1350 1351 INTEGER :: nit 1352 INTEGER :: nsd 1353 INTEGER :: nsteps 1354 ! 1355 REAL(DP) :: tolw 1356 ! 1357 LOGICAL :: adapt 1358 ! 1359 INTEGER :: calwf 1360 INTEGER :: nwf 1361 INTEGER :: wffort 1362 ! 1363 LOGICAL :: writev 1364!============================================================================== 1365!exx_wf related 1366 NAMELIST / wannier / wf_efield, wf_switch, sw_len, efx0, efy0, efz0,& 1367 efx1, efy1, efz1, wfsd, wfdt,exx_neigh,exx_poisson_eps,& 1368 exx_dis_cutoff,exx_ps_rcut_self, exx_me_rcut_self, & 1369 exx_use_cube_domain, & 1370 exx_ps_rcut_pair, exx_me_rcut_pair, vnbsp,& 1371 maxwfdt, wf_q, wf_friction, nit, nsd, nsteps, & 1372 tolw, adapt, calwf, nwf, wffort, writev 1373!=============================================================================== 1374! END manual 1375! ---------------------------------------------------------------------- 1376 1377!=----------------------------------------------------------------------------=! 1378! WANNIER_NEW Namelist Input Parameters 1379!=----------------------------------------------------------------------------=! 1380 1381 LOGICAL :: & 1382 plot_wannier = .false.,& 1383 ! if .TRUE. wannier number plot_wan_num is plotted 1384 use_energy_int = .false., & 1385 ! if .TRUE. energy interval is used to generate wannier 1386 print_wannier_coeff = .false. 1387 ! if .TRUE. 1388 INTEGER, PARAMETER :: nwanx = 50 ! max number of wannier functions 1389 INTEGER :: & 1390 nwan, &! number of wannier functions 1391 plot_wan_num = 0, &! number of wannier for plotting 1392 plot_wan_spin = 1 ! spin of wannier for plotting 1393 REAL(DP) :: & 1394 constrain_pot(nwanx,2) ! constrained potential for wannier 1395 NAMELIST / wannier_ac / plot_wannier, use_energy_int, nwan, & 1396 plot_wan_num, plot_wan_spin, constrain_pot, print_wannier_coeff 1397 1398! END manual 1399! ---------------------------------------------------------------------- 1400 1401 1402! ---------------------------------------------------------------- 1403! BEGIN manual 1404! 1405!=----------------------------------------------------------------------------=! 1406! CARDS parameters 1407!=----------------------------------------------------------------------------=! 1408! 1409! Note: See file read_cards.f90 for card syntax and usage 1410! 1411! ATOMIC_SPECIES 1412! 1413 CHARACTER(len=3) :: atom_label(nsx) = 'XX' ! label of the atomic species being read 1414 CHARACTER(len=80) :: atom_pfile(nsx) = 'YY' ! pseudopotential file name 1415 REAL(DP) :: atom_mass(nsx) = 0.0_DP ! atomic mass of the i-th atomic species 1416 ! in atomic mass units: 1 a.m.u. = 1822.9 a.u. = 1.6605 * 10^-27 kg 1417 LOGICAL :: taspc = .false. 1418 LOGICAL :: tkpoints = .false. 1419 LOGICAL :: tforces = .false. 1420 LOGICAL :: tocc = .false. 1421 LOGICAL :: tcell = .false. 1422 LOGICAL :: tionvel = .false. 1423 LOGICAL :: tconstr = .false. 1424 LOGICAL :: tksout = .false. 1425 LOGICAL :: ttemplate = .false. 1426 LOGICAL :: twannier = .false. 1427 1428! 1429! ATOMIC_POSITIONS 1430! 1431 REAL(DP), ALLOCATABLE :: rd_pos(:,:) ! unsorted positions from input 1432 INTEGER, ALLOCATABLE :: sp_pos(:) 1433 INTEGER, ALLOCATABLE :: rd_if_pos(:,:) 1434 INTEGER, ALLOCATABLE :: na_inp(:) 1435 LOGICAL :: tapos = .false. 1436 LOGICAL :: lsg = .false. 1437 CHARACTER(len=80) :: atomic_positions = 'crystal' 1438 ! atomic_positions = 'bohr' | 'angstrom' | 'crystal' | 'alat' 1439 ! select the units for the atomic positions being read from stdin 1440 1441! 1442! ION_VELOCITIES 1443! 1444 REAL(DP), ALLOCATABLE :: rd_vel(:,:) ! unsorted velocities from input 1445 INTEGER, ALLOCATABLE :: sp_vel(:) 1446 LOGICAL :: tavel = .false. 1447! 1448! ATOMIC_FORCES 1449! 1450 REAL(DP), ALLOCATABLE :: rd_for(:,:) ! external forces applied to single atoms 1451 1452! 1453! KPOINTS 1454! 1455! ... k-points inputs 1456 LOGICAL :: tk_inp = .false. 1457 REAL(DP), ALLOCATABLE :: xk(:,:), wk(:) 1458 INTEGER :: nkstot = 0, nk1 = 0, nk2 = 0, nk3 = 0, k1 = 0, k2 = 0, k3 = 0 1459 CHARACTER(len=80) :: k_points = 'gamma' 1460 ! k_points = 'automatic' | 'crystal' | 'tpiba' | 'gamma'* 1461 ! k_points = 'crystal_b' | 'tpiba_b' 1462 ! select the k points mesh 1463 ! 'automatic' k points mesh is generated automatically 1464 ! with Monkhorst-Pack algorithm 1465 ! 'crystal' k points mesh is given in stdin in scaled units 1466 ! 'tpiba' k points mesh is given in stdin in units of ( 2 PI / alat ) 1467 ! 'gamma' only gamma point is used ( default in CPMD simulation ) 1468 ! _b means that a band input is given. The weights is a integer 1469 ! number that gives the number of points between the present point 1470 ! and the next. The weight of the last point is not used. 1471! 1472! OCCUPATIONS 1473! 1474 REAL(DP), ALLOCATABLE :: f_inp(:,:) 1475 LOGICAL :: tf_inp = .false. 1476 1477! 1478! CELL_PARAMETERS 1479! 1480 REAL(DP) :: rd_ht(3,3) = 0.0_DP 1481 CHARACTER(len=80) :: cell_units = 'none' 1482 LOGICAL :: trd_ht = .false. 1483 1484! 1485! REFERENCE_CELL_PARAMETERS 1486! 1487 REAL(DP) :: rd_ref_ht(3,3) = 0.0_DP 1488 CHARACTER(len=80) :: ref_cell_units = 'alat' 1489 LOGICAL :: ref_cell = .false. 1490 1491! 1492! CONSTRAINTS 1493! 1494 INTEGER :: nc_fields = 4 ! max number of fields that is allowed to 1495 ! define a constraint 1496 1497 INTEGER :: nconstr_inp = 0 1498 REAL(DP) :: constr_tol_inp = 1.E-6_DP 1499 ! 1500 CHARACTER(len=20), ALLOCATABLE :: constr_type_inp(:) 1501 REAL(DP), ALLOCATABLE :: constr_inp(:,:) 1502 REAL(DP), ALLOCATABLE :: constr_target_inp(:) 1503 LOGICAL, ALLOCATABLE :: constr_target_set(:) 1504 1505! 1506! KOHN_SHAM 1507! 1508 INTEGER, ALLOCATABLE :: iprnks( :, : ) 1509 INTEGER :: nprnks( nspinx ) = 0 1510 ! logical mask used to specify which kohn sham orbital should be 1511 ! written to files 'KS.' 1512 1513! 1514! PLOT_WANNIER 1515! 1516 1517 INTEGER, PARAMETER :: nwf_max = 1000 1518 ! 1519 INTEGER :: wannier_index( nwf_max ) 1520 1521! 1522! WANNIER_NEW 1523! 1524 TYPE (wannier_data) :: wan_data(nwanx,2) 1525 1526 1527! END manual 1528! ---------------------------------------------------------------------- 1529 1530 LOGICAL :: xmloutput = .false. 1531 ! if .true. PW produce an xml output 1532 1533CONTAINS 1534! 1535!---------------------------------------------------------------------------- 1536SUBROUTINE reset_input_checks() 1537 !----------------------------------------------------------------------------- 1538 ! 1539 ! ... This routine sets to .false. flags used to check whether some variables 1540 ! ... have been read. If called before reading, allows to read a different 1541 ! ... input file without triggering bogus error messages - useful for NEB 1542 ! 1543 IMPLICIT NONE 1544 ! 1545 tapos = .false. 1546 tkpoints = .false. 1547 taspc = .false. 1548 twannier = .false. 1549 tconstr = .false. 1550 tforces = .false. 1551 tocc = .false. 1552 tksout = .false. 1553 tionvel = .false. 1554 tcell = .false. 1555 ! 1556 END SUBROUTINE reset_input_checks 1557 ! 1558 ! 1559 !----------------------------------------------------------------------------- 1560 SUBROUTINE allocate_input_ions( ntyp, nat ) 1561 !----------------------------------------------------------------------------- 1562 ! 1563 INTEGER, INTENT(in) :: ntyp, nat 1564 ! 1565 IF ( allocated( rd_pos ) ) DEALLOCATE( rd_pos ) 1566 IF ( allocated( sp_pos ) ) DEALLOCATE( sp_pos ) 1567 IF ( allocated( rd_if_pos ) ) DEALLOCATE( rd_if_pos ) 1568 IF ( allocated( na_inp ) ) DEALLOCATE( na_inp ) 1569 IF ( allocated( rd_vel ) ) DEALLOCATE( rd_vel ) 1570 IF ( allocated( sp_vel ) ) DEALLOCATE( sp_vel ) 1571 IF ( allocated( rd_for ) ) DEALLOCATE( rd_for ) 1572 ! 1573 ALLOCATE( rd_pos( 3, nat ) ) 1574 ALLOCATE( sp_pos( nat) ) 1575 ALLOCATE( rd_if_pos( 3, nat ) ) 1576 ALLOCATE( na_inp( ntyp) ) 1577 ALLOCATE( rd_vel( 3, nat ) ) 1578 ALLOCATE( sp_vel( nat) ) 1579 ALLOCATE( rd_for( 3, nat ) ) 1580 ! 1581 rd_pos = 0.0_DP 1582 sp_pos = 0 1583 rd_if_pos = 1 1584 na_inp = 0 1585 rd_vel = 0.0_DP 1586 sp_vel = 0 1587 rd_for = 0.0_DP 1588 ! 1589 RETURN 1590 ! 1591 END SUBROUTINE allocate_input_ions 1592 1593 !----------------------------------------------------------------------------- 1594 SUBROUTINE allocate_input_constr() 1595 !----------------------------------------------------------------------------- 1596 ! 1597 IF ( allocated( constr_type_inp ) ) DEALLOCATE( constr_type_inp ) 1598 IF ( allocated( constr_inp ) ) DEALLOCATE( constr_inp ) 1599 IF ( allocated( constr_target_inp ) ) DEALLOCATE( constr_target_inp ) 1600 IF ( allocated( constr_target_set ) ) DEALLOCATE( constr_target_set ) 1601 ! 1602 ALLOCATE( constr_type_inp( nconstr_inp ) ) 1603 ALLOCATE( constr_target_inp( nconstr_inp ) ) 1604 ALLOCATE( constr_target_set( nconstr_inp ) ) 1605 ! 1606 ALLOCATE( constr_inp( nc_fields, nconstr_inp ) ) 1607 ! 1608 constr_type_inp = ' ' 1609 constr_inp = 0.0_DP 1610 constr_target_inp = 0.0_DP 1611 constr_target_set = .false. 1612 ! 1613 RETURN 1614 ! 1615 END SUBROUTINE allocate_input_constr 1616 1617 !----------------------------------------------------------------------------- 1618 SUBROUTINE allocate_input_iprnks( nksx, nspin ) 1619 !----------------------------------------------------------------------------- 1620 ! 1621 INTEGER, INTENT(in) :: nksx, nspin 1622 ! 1623 IF( allocated( iprnks ) ) DEALLOCATE( iprnks ) 1624 ! 1625 ALLOCATE( iprnks( max( 1, nksx), nspin ) ) 1626 ! 1627 iprnks = 0 1628 ! 1629 RETURN 1630 ! 1631 END SUBROUTINE allocate_input_iprnks 1632 1633 !----------------------------------------------------------------------------- 1634 SUBROUTINE deallocate_input_parameters() 1635 !----------------------------------------------------------------------------- 1636 ! 1637 IF ( allocated( xk ) ) DEALLOCATE( xk ) 1638 IF ( allocated( wk ) ) DEALLOCATE( wk ) 1639 IF ( allocated( rd_pos ) ) DEALLOCATE( rd_pos ) 1640 IF ( allocated( sp_pos ) ) DEALLOCATE( sp_pos ) 1641 IF ( allocated( rd_if_pos ) ) DEALLOCATE( rd_if_pos ) 1642 IF ( allocated( na_inp ) ) DEALLOCATE( na_inp ) 1643 IF ( allocated( rd_vel ) ) DEALLOCATE( rd_vel ) 1644 IF ( allocated( sp_vel ) ) DEALLOCATE( sp_vel ) 1645 IF ( allocated( rd_for ) ) DEALLOCATE( rd_for ) 1646 ! 1647 IF ( allocated( f_inp ) ) DEALLOCATE( f_inp ) 1648 ! 1649 IF ( allocated( constr_type_inp ) ) DEALLOCATE( constr_type_inp ) 1650 IF ( allocated( constr_inp ) ) DEALLOCATE( constr_inp ) 1651 IF ( allocated( constr_target_inp ) ) DEALLOCATE( constr_target_inp ) 1652 IF ( allocated( constr_target_set ) ) DEALLOCATE( constr_target_set ) 1653 ! 1654 IF ( allocated( iprnks ) ) DEALLOCATE( iprnks ) 1655 ! 1656 RETURN 1657 ! 1658 END SUBROUTINE deallocate_input_parameters 1659 ! 1660!=----------------------------------------------------------------------------=! 1661! 1662END MODULE input_parameters 1663! 1664!=----------------------------------------------------------------------------=! 1665