Input data format: { } = optional, [ ] = it depends All quantities whose dimensions are not explicitly specified are in RYDBERG ATOMIC UNITS =============================================================================== &CONTROL ... / &SYSTEM ... / &ELECTRONS ... / [ &IONS ... / ] [ &CELL ... / ] [ &PHONON ... / ] ATOMIC_SPECIES X Mass_X PseudoPot_X Y Mass_Y PseudoPot_Y Z Mass_Z PseudoPot_Z ATOMIC_POSITIONS { alat | bohr | crystal | angstrom } in all cases except calculation = 'neb' or 'smd' : X 0.0 0.0 0.0 {if_pos(1) if_pos(2) if_pos(3)} Y 0.5 0.0 0.0 Z O.0 0.2 0.2 if calculation = 'neb' .OR. 'smd' : first_image X 0.0 0.0 0.0 {if_pos(1) if_pos(2) if_pos(3)} Y 0.5 0.0 0.0 Z O.0 0.2 0.2 { intermediate_image 1 X 0.0 0.0 0.0 Y 0.9 0.0 0.0 Z O.0 0.2 0.2 intermediate_image ... X 0.0 0.0 0.0 Y 0.9 0.0 0.0 Z O.0 0.2 0.2 } last_image X 0.0 0.0 0.0 Y 0.7 0.0 0.0 Z O.0 0.5 0.2 K_POINTS { tpiba | automatic | crystal | gamma } if (gamma) nothing to read if (automatic) nk1, nk2, nk3, k1, k2, k3 if (not automatic) nks xk_x, xk_y, xk_z, wk [ CELL_PARAMETERS { cubic | hexagonal } a(1,1) a(2,1) a(3,1) a(1,2) a(2,2) a(3,2) a(1,3) a(2,3) a(3,3) ] [ OCCUPATIONS f_inp(1,1) f_inp(2,1) f_inp(3,1) ... f_inp(10,1) f_inp(11,1) f_inp(12,1) ... f_inp(nbnd,1) [ f_inp(1,2) f_inp(2,2) f_inp(3,2) ... f_inp(10,2) f_inp(11,2) f_inp(12,2) ... f_inp(nbnd,2) ] ] [ CLIMBING_IMAGES list of images, separated by a comma ] [ CONSTRAINTS nconstr constr_tol constr_type(.) constr(1,.) constr(2,.) { constr_target(.) } ] [ COLLECTIVE_VARS nconstr constr_tol constr_type(.) constr(1,.) constr(2,.) { constr_target(.) } ] =============================================================================== NAMELIST &CONTROL calculation CHARACTER a string describing the task to be performed: 'scf', 'nscf', 'bands', 'phonon', 'relax', 'md', 'vc-relax', 'vc-md', 'neb', 'smd', 'metadyn' (vc = variable-cell). Default: 'scf' title CHARACTER reprinted on output. Default: ' ' verbosity CHARACTER 'high' | 'default' | 'low' | 'minimal' restart_mode CHARACTER 'from_scratch' : from scratch ( default ) NEB and SMD only: the starting path is obtained with a linear interpolation between the images specified in the ATOMIC_POSITIONS card. Note that in the linear interpolation periodic boundary conditions ARE NON USED. 'restart' : from previous interrupted run wf_collect LOGICAL ( default = .FALSE. ) This flag controls the way wavefunctions are stored to disk : .TRUE. collect wavefunctions from all processors and store them into the output data directory outdir/prefix.save .FALSE. do not collect wavefunctions, leave them in temporary local files (one per processor). The resulting format will be readable only by jobs running on the same number of processors and pools. Useful if you do not need the wavefunction or if you want to reduce the I/O or the disk occupancy. nstep INTEGER number of ionic + electronic steps default: 1 if calculation = 'scf', 'nscf', 'bands' 0 if calculation = 'neb', 'smd' 50 for the other cases iprint INTEGER band energies are written every iprint iterations default: write only at convergence tstress LOGICAL calculate stress. Set to .TRUE. if calculation='vc-md' tprnfor LOGICAL print forces. Set to .TRUE. if calculation='relax','md','vc-md' dt REAL ( default = 20.D0 ) time step for molecular dynamics, in Rydberg atomic units (1 a.u.=4.8378 * 10^-17 s : beware, CP and FPMD codes use Hartree atomic units, half that much!!!) outdir CHARACTER ( default = value of the ESPRESSO_TMPDIR environment variable if set; current directory ('./') otherwise ) input, temporary, output files are found in this directory, see also 'wfcdir' wfcdir CHARACTER ( by default same as outdir ) this directory specifies where to store files generated by each processor (*.wfc{N}, *.igk{N}, etc.). The idea here is to be able to separately store the largest files, while the files necessary for restarting still go into 'outdir' (for now only works for stand alone PW ) prefix CHARACTER ( default = 'pwscf' ) prepended to input/output filenames: prefix.wfc, prefix.rho, etc. lkpoint_dir LOGICAL (default=.true.) If .false. it does not open a subdirectory for each k_point in the prefix.save directory. max_seconds REAL ( default : 1.D+7, or 150 days, i.e. no time limit ) jobs stops after max_seconds CPU time etot_conv_thr REAL ( default = 1.0D-4 ) convergence threshold on total energy (a.u) for ionic minimization: the convergence criterion is satisfied when the total energy changes less than etot_conv_thr between two consecutive scf steps. See also forc_conv_thr - both criteria must be satisfied forc_conv_thr REAL ( default = 1.0D-3 ) convergence threshold on forces (a.u) for ionic minimization: the convergence criterion is satisfied when all components of all forces are smaller than forc_conv_thr. See also etot_conv_thr - both criteria must be satisfied disk_io CHARACTER ( default='default') Specifies the amount of disk I/O activity 'high': save all data at each SCF step 'default': save wavefunctions at each SCF step unless there is a single k-point per process 'low' : store wfc in memory, save only at the end 'none': do not save wfc, not even at the end If restarting from an interrupted calculation, the code will try to figure out what is available on disk. The more you write, the more complete the restart will be. pseudo_dir CHARACTER ( default = value of the $ESPRESSO_PSEUDO environment variable if set; '$HOME/espresso/pseudo/' otherwise ) directory containing pseudopotential files tefield LOGICAL ( default = .FALSE. ) If .TRUE. a sawlike potential simulating an electric field is added to the bare ionic potential. See variables edir, eamp, emaxprog, eopreg for the form and size of the added potential. dipfield LOGICAL ( default = .FALSE. ) If .TRUE. and tefield=.TRUE. a dipole correction is also added to the bare ionic potential - implements the recipe of L. Bengtsson, PRB 59, 12303 (1999). See variables edir, emaxprog, eopreg for the form of the correction, that must be used only in a slab geometry, for surface calculations, with the discontinuity in the empty space lelfield LOGICAL ( default = .FALSE. ) If .TRUE. a homogeneous finite electric field described through the modern theory of the polarization is applied. This is different from "tefield=.true." ! lberry LOGICAL (default = .FALSE.) If .TRUE. perform a Berry phase calculation See the header of PW/bp_c_phase.f90 for documentation gdir INTEGER For Berry phase calculation: direction of the k-point strings in reciprocal space. Allowed values: 1, 2, 3 1=first, 2=second, 3=third reciprocal lattice vector For calculations with finite electric fields (lelfield==.true.), gdir is the direction of the field nppstr INTEGER For Berry phase calculation: number of k-points to be calculated along each symmetry-reduced string The same for calculation with finite electric fields (lelfield==.true.) nberrycyc INTEGER ( default = 1 ) In the case of a finite electric field ( lelfield == .TRUE. ) it defines the number of iterations for converging the wavefunctions in the electric field Hamiltonian, for each external iteration on the charge density =============================================================================== NAMELIST &SYSTEM ibrav INTEGER bravais-lattice index - see at the end of this file celldm(i) REAL, DIMENSION(6) crystallographic constants - see at the end of this file alat = celldm(1) is the lattice parameter "a" (in BOHR) only needed celldm (depending on ibrav) must be specified If ibrav=0 only alat = celldm(1) is used (if present) a, b, c, cosab, cosac, cosbc: REAL traditional crystallographic constants (a,b,c in ANGSTROM, cosab = cosine of the angle between axis a and b specify either these OR celldm but NOT both If ibrav=0 only alat = a is used (if present) nat INTEGER number of atoms in the unit cell - must be specified ntyp INTEGER number of types of atoms in the unit cell - must be specified nbnd INTEGER number of electronic states (bands) to be calculated. Default: for an insulator, nbnd = (number of valence bands) (nbnd=nelec/2, see below for nelec) for a metal, 20% more (minimum 4 more) Note that in spin-polarized calculations the number of k-point, not the number of bands per k-point, is doubled nelec REAL number of electron in the unit cell (may be noninteger if you wish) Default: the same as ionic charge (neutral cell) A compensating jellium background is inserted to remove divergences if the cell is not neutral tot_charge INTEGER ( default = 0 ) total system charge. Used only if nelec is unspecified, otherwise it is ignored. ecutwfc REAL kinetic energy cutoff (Ry) for wavefunctions (must be specified) ecutrho REAL ( default = 4 * ecutwfc ) kinetic energy cutoff (Ry) for charge density and potential May be larger ( for ultrasoft PP ) or somewhat smaller ( but not much smaller ) than the default value. Note that if you have norm-conserving PP only, setting it to a larger value than the default is a waste of time. nr1,nr2,nr3 INTEGER three-dimensional FFT mesh (hard grid) for charge density (and scf potential). If not specified the grid is calculated based on the cutoff for charge density (see also "ecutrho") nr1s,nr2s,nr3s INTEGER three-dimensional mesh for wavefunction FFT and for the smooth part of charge density ( smooth grid ). Coincides with nr1, nr2, nr3 if ecutrho = 4 * ecutwfc ( default ) nosym LOGICAL ( default = .FALSE. ) if (.TRUE.) symmetry is not used. Note that a k-point grid provided in input is used "as is"; an automatically generated k-point grid will contain only points in the irreducible BZ of the lattice. Use with care in low-symmetry large cells if you cannot afford a k-point grid with the correct symmetry. occupations CHARACTER 'smearing': gaussian smearing for metals requires a value for degauss 'tetrahedra' : for metals and DOS calculation (see PRB49, 16223 (1994)) Requires uniform grid of k-points, automatically generated (see below) Not suitable (because not variational) for force/optimization/dynamics calculations 'fixed' : for insulators with a gap 'from_input' : The occupation are read from input file. Presently works only with one k-point (LSDA allowed). degauss REAL ( default = 0.D0 Ry ) value of the gaussian spreading (Ry) for brillouin-zone integration in metals. smearing CHARACTER 'gaussian', 'gauss': ordinary Gaussian spreading (Default) 'methfessel-paxton', 'm-p', 'mp': Methfessel-Paxton first-order spreading (see PRB 40, 3616 (1989)). 'marzari-vanderbilt', 'cold', 'm-v', 'mv': Marzari-Vanderbilt cold smearing (see PRL 82, 3296 (1999)) 'fermi-dirac', 'f-d', 'fd': smearing with Fermi-Dirac function nspin INTEGER nspin = 1 : non-polarized calculation (default) nspin = 2 : spin-polarized calculation, LSDA (magnetization along z axis) nspin = 4 : spin-polarized calculation, noncollinear (magnetization in generic direction) DO NOT specify nspin in this case; specify "noncolin=.TRUE." instead noncolin LOGICAL if .true. the program will perform a noncollinear calculation. DEFAULT: .false. starting_magnetization(i) REAL starting spin polarization (values between -1 and 1) on atomic type 'i' in a spin-polarized calculation. Breaks the symmetry and provides a starting point for self-consistency. The default value is zero, BUT a value MUST be specified for AT LEAST one atomic type in spin polarized calculations. Note that if start from zero initial magnetization, you will get zero final magnetization in any case. If you desire to start from an antiferromagnetic state, you may need to define two different atomic species corresponding to sublattices of the same atomic type. If you fix the magnetization with "nelup/neldw" or with "multiplicity" or with "tot_magnetization", you should not specify starting_magnetization. If you are restarting from a previous run, or from an interrupted run, starting_magnetization is ignored. nelup, neldw REAL number of spin-up and spin-down electrons, respectively Note that this fixes the final value of the magnetization. The sum must yield nelec that must also be specified explicitly in this case. Not valid for spin-unpolarized or noncollinear calculations, only for LSDA. Obsolescent: use multiplicity or tot_magnetization instead. multiplicity INTEGER ( default = 0 [unspecified] ) spin multiplicity (2s+1). 1 is singlet, 2 for doublet etc. Note that this fixes the final value of the magnetization. if unspecified or a non-zero value is specified in nelup/neldw then multiplicity variable is ignored. Do not specify both multiplicity and tot_magnetization. tot_magnetization INTEGER ( default = -1 [unspecified] ) majority spin - minority spin (nelup - neldw). if unspecified or a non-zero value is specified in nelup/neldw then tot_magnetization variable is ignored. Do not specify both multiplicity and tot_magnetization. YES, there is redundancy! nelup/neldw are enough to specify the spin state. However these variables are not very convenient and will be eliminated from the input in future versions. It is recommended to use either 'multiplicity' or equivalently 'tot_magnetization' to specify the spin state. ecfixed REAL ( default = 0.0 ) qcutz REAL ( default = 0.0 ) q2sigma REAL ( default = 0.1 ) parameters for modified functional to be used in variable-cell molecular dynamics (or in stress calculation). "ecfixed" is the value (in Rydberg) of the constant-cutoff; "qcutz" and "q2sigma" are the height and the width (in Rydberg) of the energy step for reciprocal vectors whose square modulus is greater than "ecfixed". In the kinetic energy, G^2 is replaced by G^2 + qcutz * (1 + erf ( (G^2 - ecfixed)/q2sigma) ) See: M. Bernasconi et al, J. Phys. Chem. Solids 56, 501 (1995) xc_type CHARACTER Exchange-correlation functional Presently unused: XC functional is read from PP files lda_plus_u LOGICAL ( default = .FALSE.) Hubbard_U(I) REAL ( default = 0.D0 for all species) Hubbard_alpha(I) REAL ( default = 0.D0 for all species) parameters for LDA+U calculations If lda_plus_u = .TRUE. you must specify, for species I, the parameters U and (optionally) alpha of the Hubbard model (both in eV). See: Anisimov, Zaanen, and Andersen, PRB 44, 943 (1991); Anisimov et al., PRB 48, 16929 (1993); Liechtenstein, Anisimov, and Zaanen, PRB 52, R5467 (1994); Cococcioni and de Gironcoli, PRB 71, 035105 (2005). IMPORTANT: LDA+U works only for a few selected elements. Modify PW/set_hubbard_l.f90 and PW/tabd.f90 if you plan to use LDA+U with an element that is not configured there. starting_ns_eigenvalue(m,ispin,I) REAL (default = -1.d0 that means NOT SET) In the first iteration of an LDA+U run it overwrites the m-th eigenvalue of the ns occupation matrix for the ispin component of atomic species I. Leave unchanged eigenvalues that are not set. This is useful to suggest the desired orbital occupations when the default choice takes another path. U_projection_type CHARACTER (default='atomic') Only active when lda_plus_U is .true., specifies the type of projector on localized orbital to be used in the LDA+U scheme. Currently available choices: 'atomic': use atomic wfc's (as they are) to build the projector 'ortho-atomic': use Lowdin orthogonalized atomic wfc's 'norm-atomic': Lowdin normalization of atomic wfc. Keep in mind: atomic wfc are not orthogonalized in this case. This is a "quick and dirty" trick to be used when atomic wfc from the pseudopotential are not normalized (and thus produce occupation whose value exceeds unity). If orthogonalized wfc are not needed always try 'atomic' first. 'file': use the information from file "prefix".atwfc that must have been generated previously, for instance by pmw.x (see PP/poormanwannier.f90 for details) NB: forces and stress currently implemented only for the 'atomic' choice. edir INTEGER The direction of the electric field or dipole correction is parallel to the bg(:,edir) reciprocal lattice vector, so the potential is constant in planes defined by FFT grid points; edir = 1, 2 or 3. Used only if tefield is .TRUE. emaxpos REAL ( default = 0.5D0 ) Position of the maximum of the sawlike potential along crystal axis "edir", within the unit cell (see below), 0 < emaxpos < 1 Used only if tefield is .TRUE. eopreg REAL( default = 0.1D0 ) Zone in the unit cell where the sawlike potential decreases. ( see below, 0 < eopreg < 1 ). Used only if tefield is .TRUE. eamp REAL ( default = 0.001 a.u. ) Amplitude of the electric field (in a.u. = 51.44 10^10 V/m ) The sawlike potential increases with slope "eamp" in the region from (emaxpos+eopreg-1) to (emaxpos), then decreases to 0 until (emaxpos+eopreg), in units of the crystal vector "edir". Used only if tefield is .TRUE. angle1(i) REAL The angle expressed in degrees between the initial magnetization and the z-axis. For noncollinear calculations only. I runs over the atom types. angle2(i) REAL The angle expressed in degrees between the projection of the initial magnetization on x-y plane and the x-axis. For noncollinear calculations only. constrained_magnetization CHARACTER ( defalt = 'none' ) Used to perform constrained calculations in magnetic systems Currently available choices: `none` : no constraint `total`: total magnetization is constrained If nspin=4 (noncolin=.True.) constraint is imposed by adding a penalty functional to the total energy: - LAMBDA * SUM_{i} ( magnetization(i) - fixed_magnetization(i) )**2 where the sum over i runs over the three components of the magnetization. Lambda is a real number (see below). If nspin=2 constraint is imposed by defining two Fermi energies for spin up and down. Only fixed_magnetization(3) can be defined in this case. `atomic`: atomic magnetization are constrained to the defined starting magnetization adding a penalty - LAMBDA * SUM_{i,itype} ( magnetic_moment(i,itype) - mcons(i,itype) )**2 where i runs over the components (1-3) and itype over the types (1-ntype). mcons(:,:) array is defined from starting_magnetization, angle1 and angle2 variables. lambda is a real number `total direction`: the angle theta of the total magnetization with the z axis (theta = fixed_magnetization(3)) is constrained: - LAMBDA * ( magnetization(1) - magnetization(3)*tan(theta) )**2 `atomic direction`: not all the components of the atomic magnetic moment are constrained but only the cosine of angle1, and the penalty functional is: - LAMBDA * SUM_{itype} ( mag_mom(3,itype)/mag_mom_tot - cos(angle1(ityp) )**2 fixed_magnetization(3) REAL (default 0.d0) value of the total magnetization to be maintained fixed when constrained_magnetization='total' B_field(3) REAL (default = 0.d0) A fixed magnetic field defined by the vector B_field is added to the exchange and correlation magnetic field. The three components of the magnetic field are given in Ry. Only B_field(3) can be used if nspin=2. In all calculations with a finite magnetic field, we print the total energy WITHOUT the B dot M term. In the calculations with the penalty functional we write only the total energy, NOT the penalty functional. lambda REAL parameter used for constrained_magnetization calculations NB: LAMBDA is reduced in the first iterations and is increased slowly up to the input value. report INTEGER It's the number of iterations after which the program write all the atomic magnetic moments. lspinorb LOGICAL if .TRUE. the noncollinear code can use a pseudopotential with spin-orbit. assume_isolated LOGICAL ( default = .FALSE. ) if .TRUE. the system is assumed to be isolated (a molecule or cluster in a supercell) and the Makov-Payne correction to the total energy is computed. An estimate of the vacuum level is also calculated so that eigenvalues can be properly aligned. =============================================================================== NAMELIST &ELECTRONS electron_maxstep INTEGER ( default = 100 ) maximum number of iterations in a scf step conv_thr REAL ( default = 1.D-6 ) Convergence threshold for selfconsistency: estimated energy error < conv_thr mixing_mode CHARACTER 'plain' : charge density Broyden mixing ( default ) 'TF' : as above, with simple Thomas-Fermi screening (for highly homogeneous systems) 'local-TF': as above, with local-density-dependent TF screening (for highly inhomogeneous systems) mixing_beta REAL ( default = 0.7D0 ) mixing factor for self-consistency mixing_ndim INTEGER ( default = 8) number of iterations used in mixing scheme mixing_fixed_ns INTEGER ( default = 0 ) For LDA+U : number of iterations with fixed ns ( ns is the atomic density appearing in the Hubbard term ) diagonalization CHARACTER 'david': Davidson iterative diagonalization with overlap matrix (default) 'cg' : conjugate-gradient-like band-by-band diagonalization 'diis' : DIIS-like diagonalization - PRESENTLY DISABLED 'david+para': like david, but use parallel diagonalization, after testing if it is convenient (for speed), iteration loop use replicated data. 'david+distpara': like david but use fully distributed memory iteration loop. Allocated memory scale down with the number of procs. Procs involved in diagonalization can be changed with input parameter "ortho_para". On multi core/CPUs often it is convenient to let only one core per CPU to work on liner algebra. ortho_para INTEGER ( default = 0 ) meaningful only if diagonalization = 'david+para' or diagonalization = 'david+distpara' and with parallel executable. The number of processors to be used for parallel diagonalization. With the default value ( 0 ) the code try to use all the processors. Actually the number of processors used is the largest square number less or equal to ortho_para. diago_thr_init REAL Convergence threshold for the first iterative diagonalization (the check is on eigenvalue convergence). For scf calculations, the default is 1.D-2 if starting from a superposition of atomic orbitals; 1.D-5 if starting from a charge density. During self consistency the threshold (ethr) is automatically reduced when approaching convergence. For non-scf calculations, this is the threshold used in the iterative diagonalization. The default is conv_thr / nelec. For 'phonon' calculations, diago_thr_init is ignored: the threshold is always set to conv_thr / nelec . diago_cg_maxiter INTEGER For conjugate gradient diagonalization: max number of iterations diago_david_ndim INTEGER ( default = 4 ) For Davidson diagonalization: dimension of workspace (number of wavefunction packets, at least 2 needed). A larger value may yield a faster algorithm but uses more memory diago_diis_ndim INTEGER ( default = 3 ) For DIIS: dimension of the reduced space. diago_full_acc LOGICAL ( default = .FALSE. ) If .TRUE. all the empty states are diagonalized at the same level of accuracy of the occupied ones. Otherwise the empty states are diagonalized using a larger threshold (this should not affect total energy, forces, and other ground-state properties). efield REAL ( default = 0.D0 ) For finite electric field calculations (lelfield == .TRUE.), it defines the intensity of the field in a.u. startingpot CHARACTER 'atomic': starting potential from atomic charge superposition ( default for scf, *relax, *md, neb, smd ) 'file' : start from existing "prefix".pot file ( default, only possibility for nscf, bands, phonon ) startingwfc CHARACTER 'atomic': start from superposition of atomic orbitals ( default ) If not enough atomic orbitals are available, fill with random numbers the remaining wfcs 'random': start from random wfcs 'file': start from a wavefunction file tqr LOGICAL (default=.FALSE.) If .true., use the (VERY EXPERIMENTAL) real-space algorithm for augmentation charges in ultrasoft pseudopotentials. Must faster execution of ultrasoft-related calculations, but numerically less accurate than the default algorithm. Use with care and after testing! =============================================================================== NAMELIST &IONS ( only if calculation = 'relax', 'md', 'vc-relax', 'vc-md', 'neb', 'smd' ) ion_dynamics CHARACTER specify the type of ionic dynamics. For constrained dynamics or constrained optimisations add the CONSTRAINTS card (when the card is present the SHAKE algorithm is automatically used). For different type of calculation different possibilities are allowed and different default values apply: CASE ( calculation = 'relax' ) 'bfgs' : (default) a new BFGS quasi-newton algorithm, based on the trust radius procedure, is used for structural relaxation (experimental) 'damp' : use damped (quick-min Verlet) dynamics for structural relaxation CASE ( calculation = 'md' ) 'verlet' : (default) use Verlet algorithm to integrate Newton's equation 'langevin' ion dynamics is over-damped Langevin CASE ( calculation = 'vc-relax' ) 'damp' : (default) use damped (Beeman) dynamics for structural relaxation CASE ( calculation = 'vc-md' ) 'beeman' : (default) use Beeman algorithm to integrate Newton's equation phase_space CHARACTER ( defauld = 'full' ) 'full' : the full phase-space is used for the ionic dynamics. 'coarse-grained' : a coarse-grained phase-space, defined by a set of constraints, is used for the ionic dynamics (used for calculation of free-energy barriers) pot_extrapolation CHARACTER ( default = "atomic" ) used to extrapolate the potential from preceding ionic steps 'none' : no extrapolation 'atomic' : extrapolate the potential as if it was a sum of atomic-like orbitals 'first_order' : extrapolate the potential with first-order formula 'second_order': as above, with second order formula wfc_extrapolation CHARACTER ( default = "none" ) used to extrapolate the wavefunctions from preceding ionic steps 'none' : no extrapolation 'first_order' : extrapolate the wave-functions with first-order formula - NOT IMPLEMENTED WITH USPP 'second_order': as above, with second order formula NOT IMPLEMENTED WITH USPP remove_rigid_rot LOGICAL ( default = .FALSE. ) this keyword is useful when simulating the dynamics and/or the thermodynamics of an isolated system. If set to true the total torque of the internal forces is set to zero by adding new forces that compensate the spurious interaction with the periodic images. This allowes for the use of smaller supercells. BEWARE: since the potential energy is no longer consistent with the forces (it still contains the spurious interaction with the repeated images), the total energy is not conserved anymore. However the dynamical and thermodynamical properties should be in closer agreement with those of an isolated system. Also the final energy of a structural relaxation will be higher, but the relaxation itself should be faster. ! ! ... keywords used for molecular dynamics ! ion_temperature CHARACTER 'rescaling' control ionic temperature via velocity rescaling (first method) see parameters "tempw" and "tolp" This is the only method implemented in VC-MD 'rescale-v' control ionic temperature via velocity rescaling (second method) see parameters "tempw" and "nraise" 'rescale-T' control ionic temperature via velocity rescaling (third method) see parameter "delta_t" 'reduce-T' reduce ionic temperature every "nraise" steps by the (negative) value "delta_t" 'berendsen' control ionic temperature using "soft" velocity rescaling - see parameters "tempw" and "nraise" 'andersen' control ionic temperature using Andersen thermostat see parameters "tempw" and "nraise" 'not_controlled' (default) ionic temperature is not controlled tempw REAL ( default = 300.D0 ) starting temperature (Kelvin) in MD runs target temperature for most thermostats tolp REAL ( default = 100.D0 ) tolerance for velocity rescaling. Velocities are rescaled if the run-averaged and target temperature differ more than tolp. delta_t REAL ( default = 1.D0 ) if ion_temperature='rescale-T': at each step the instantaneous temperature is multiplied by delta_t; this is done rescaling all the velocities. if ion_temperature='reduce-T': every 'nraise' steps the instantaneous temperature is reduced by -delta_T (.e. delta_t is added to the temperature) The instantaneous temperature is calculated at the end of every ionic move and BEFORE rescaling. This is the temperature reported in the main output. For delta_t < 0, the actual average rate of heating or cooling should be roughly C*delta_t/(nraise*dt) (C=1 for an ideal gas, C=0.5 for a harmonic solid, theorem of energy equipartition between all quadratic degrees of freedom). nraise INTEGER ( default = 1 ) if ion_temperature='reduce-T': every 'nraise' steps the instantaneous temperature is reduced by -delta_T (.e. delta_t is added to the temperature) if ion_temperature='rescale-v': every 'nraise' steps the average temperature, computed from the last nraise steps, is rescaled to tempw if ion_temperature='berendsen': the "rise time" parameter is given in units of the time step: tau = nraise*dt, so dt/tau = 1/nraise if ion_temperature='andersen': the "collision frequency" parameter is given as nu=1/tau defined above, so nu*dt = 1/nraise refold_pos LOGICAL ( default = .FALSE. ) this keyword applies only in the case of molecular dynamics or damped dynamics. If true the ions are refolded at each step into the supercell. ! ! ... keywords used only in BFGS calculations ! upscale REAL ( default = 10.D0 ) max reduction factor for conv_thr during structural optimization conv_thr is automatically reduced when the relaxation approaches convergence so that forces are still accurate, but conv_thr will not be reduced to less that conv_thr / upscale bfgs_ndim INTEGER ( default = 1 ) number of old forces and displacements vectors used in the PULAY mixing of the residual vectors obtained on the basis of the inverse hessian matrix given by the BFGS algorithm. When bfgs_ndim = 1, the standard quasi-Newton BFGS method is used. (bfgs only) trust_radius_max REAL ( default = 0.8D0 ) maximum ionic displacement in the structural relaxation (bfgs only) trust_radius_min REAL ( default = 1.D-3 ) minimum ionic displacement in the structural relaxation BFGS is reset when trust_radius < trust_radius_min (bfgs only) trust_radius_ini REAL ( default = 0.5D0 ) initial ionic displacement in the structural relaxation (bfgs only) w_1, w_2 REAL ( w_1 = 0.01D0, w_2 = 0.5D0 ) parameters used in line search based on the Wolfe conditions (bfgs only) ! ! ... keywords used only in NEB and SMD calculations ! num_of_images INTEGER ( default = 0 ) number of points used to discretize the path (it must be larger than 3) opt_scheme CHARACTER ( default = "quick-min" ) specify the type of optimization scheme "sd" : steepest descent "broyden" : quasi-Newton Broyden's second method (suggested) "quick-min" : an optimisation algorithm based on the projected velovity Verlet scheme "langevin" : finite temperature langevin dynamics of the string (smd only). It is used to compute the average path and the free-energy profile. CI_scheme CHARACTER. ( default = "no-CI" ) specify the type of Climbing Image scheme "no-CI" : climbing image is not used "auto" : original CI scheme. The image highest in energy does not feel the effect of springs and is allowed to climb along the path "manual" : images that have to climb are manually selected. See also CLIMBING_IMAGES card first_last_opt LOGICAL ( default = .FALSE. ) also the first and the last configurations are optimised "on the fly" (these images do not feel the effect of the springs) temp_req REAL ( default = 0.D0 Kelvin ) temperature used for the langevin dynamics of the string. ds REAL ( default = 1.D0 ) optimisation step length ( Hartree atomic units ). If opt_scheme="broyden", ds is used as a guess for the diagonal part of the Jacobian matrix. k_max, k_min REAL ( default = 0.1D0 Hartree atomic units ) set them to use a Variable Elastic Constants scheme elastic constants are in the range [ k_min, k_max ] this is useful to rise the resolution around the saddle point path_thr REAL ( default = 0.05D0 eV / Angstrom ) the simulation stops when the error ( the norm of the force orthogonal to the path in eV/A ) is less than path_thr. use_masses LOGICAL ( default = .FALSE. ) If. TRUE. the optimisation of the path is performed using mass-weighted coordinates. use_freezing LOGICAL ( default = .FALSE. ) If. TRUE. the images are optimised according to their error: only those images with an error larger than half of the largest are optimised. The other images are kept frozen. ! ! ... keywords used only in meta-dynamics calculations ( see also the card ! ... COLLECTIVE_VARS ) ! fe_step(i) REAL ( default = 0.04 ) meta-dynamics step length (in principle different for each collective variable), defined using the same units used to define the collective variables themselves. The step also defines the spread of the Gaussian-like bias potential. g_amplitude REAL ( default = 0.005 Hartree ) Amplitude of the gaussians used in meta-dynamics. fe_nstep INTEGER ( default = 100 ) Maximum number of steps used to evaluate the potential of mean force. sw_nstep INTEGER ( default = 10 ) Number of steps used to switch to the new values of the collective variables. =============================================================================== NAMELIST &CELL ( only if calculation = 'vc-relax', 'vc-md' ) cell_dynamics CHARACTER specify the type of dynamics for the cell. For different type of calculation different possibilities are allowed and different default values apply: CASE ( calculation = 'vc-relax' ) 'none': default 'sd': steepest descent ( not implemented ) 'damp-pr': damped (Beeman) dynamics of the Parrinello-Rahman extended lagrangian 'damp-w': damped (Beeman) dynamics of the new Wentzcovitch extended lagrangian CASE ( calculation = 'vc-md' ) 'none': default 'pr': (Beeman) molecular dynamics of the Parrinello-Rahman extended lagrangian 'w': (Beeman) molecular dynamics of the new Wentzcovitch extended lagrangian press REAL ( default = 0.D0 ) target pressure [KBar] in a variable-cell md or relaxation run. wmass REAL ( default=0.75*Tot_Mass/pi**2 for Parrinello-Rahman MD, =0.75*Tot_Mass/pi**2/Omega**(2/3) for Wentzcovitch MD.) fictitious cell mass [amu] for variable-cell simulations (both 'vc-md' and 'vc-relax') cell_factor REAL ( default = 1.2D0 ) used in the construction of the pseudopotential tables. It should exceed the maximum linear contraction of the cell during a simulation. press_cov_thr REAL ( default = 0.5D0 Kbar ) convergence threshold on the pressure for variable cell relaxation ('vc-relax' : note that the other convergence thresholds for ionic relaxation apply as well). cell_dofree CHARACTER ( default = 'all' ) select which of the cell parameters should be moved all = all axis and angles are propagated volume = the cell is simply rescaled, without changing the shape x = only the x axis is moved y = only the y axis is moved z = only the z axis is moved xy = only the x and y axis are moved, angles are unchanged xz = only the x and z axis are moved, angles are unchanged yz = only the y and z axis are moved, angles are unchanged xyz = x, y and z axis are moved, angles are unchanged xyt = x1, x2, y2 (i.e. lower xy triangle of the 2 vectors) xys = x1, y1, x2, y2 (i.e. xy square of the 2 vectors) xyzt = x1, x2, y2, x3, y3, z3 (i.e. lower xyz triangle of the 3 vectors) =============================================================================== &PHONON ( only in calculation = 'phonon' ) modenum INTEGER ( default = 0 ) for single-mode phonon calculation : modenum is the index of the irreducible representation (irrep) into which the reducible representation formed by the 3*nat atomic displacements are decomposed in order to perform the phonon calculation xqq(3) REAL q-point (units 2pi/a) for phonon calculation =============================================================================== CARDS: { } = optional ------------------------------------------------------------------------------- ATOMIC_SPECIES Syntax: ATOMIC_SPECIES X(1) Mass_X(1) PseudoPot_X(ntyp) X(2) Mass_X(2) PseudoPot_X(ntyp) ... X(ntyp) Mass_X(ntyp) PseudoPot_X(ntyp) Description: X CHARACTER : label of the atom Mass_X REAL : mass of the atomic species [amu: mass of C = 12] not used if calculation='scf','nscf', 'bands', 'phonon' PseudoPot_X CHARACTER: file containing PP for this species The pseudopotential file is assumed to be in the new UPF format. If it doesn't work, the pseudopotential format is determined by the file name: *.vdb or *.van Vanderbilt US pseudopotential code *.RRKJ3 Andrea Dal Corso's code (old format) none of the above old PWscf norm-conserving format ------------------------------------------------------------------------------- ATOMIC_POSITIONS { alat | bohr | angstrom | crystal } alat : atomic positions are in cartesian coordinates, in units of the lattice parameter "a" (default) bohr : atomic positions are in cartesian coordinate, in atomic units (i.e. Bohr) angstrom: atomic positions are in cartesian coordinates, in Angstrom crystal : atomic positions are in crystal coordinates, i.e. in relative coordinates of the primitive lattice vectors (see below) - in all cases EXCEPT calculation = 'neb' or 'smd' : There are "nat" cards like the following X x y z {if_pos(1) if_pos(2) if_pos(3)} where : X Character: label of the atom as specified in ATOMIC_SPECIES x, y, z Real: atomic positions if_pos: Integer, optional ( default = 1 ): component i of the force for this atom is multiplied by if_pos(i), which must be either 0 or 1. Used to keep selected atoms and/or selected components fixed in meta-dynamics, neb, smd, MD dynamics or structural optimization run. - if calculation = 'neb' .OR. 'smd' There are at least two groups of cards, each group composed by an identifier followed by "nat" cards as specified above: identifier X x y z {if_pos(1) if_pos(2) if_pos(3)} The first group ( identifier="first_image" ) contains the first image, the last group ( identifier="last_image" ) contains the last image. There is also the possibility of specifying intermediate images; in this case their coordinates must be set between the first_image and the last_image. ( identifier="intermediate_image", followed by "nat" position cards ). Image configurations must be specified in the following order: first_image <= mandatory X 0.0 0.0 0.0 { if_pos(1) if_pos(2) if_pos(3) } Y 0.5 0.0 0.0 { if_pos(1) if_pos(2) if_pos(3) } Z 0.0 0.2 0.2 { if_pos(1) if_pos(2) if_pos(3) } intermediate_image 1 <= optional X 0.0 0.0 0.0 Y 0.9 0.0 0.0 Z 0.0 0.2 0.2 intermediate_image ... <= optional X 0.0 0.0 0.0 Y 0.9 0.0 0.0 Z 0.0 0.2 0.2 last_image <= mandatory X 0.0 0.0 0.0 Y 0.7 0.0 0.0 Z 0.0 0.5 0.2 IMPORTANT: the total number of configurations specified in the input file must be less than num_of_images (as specified in &IONS). The initial path is obtained interpolating between the specified configurations so that all images are equispaced (only the coordinates of the first and last images are not changed). ------------------------------------------------------------------------------- K_POINTS { tpiba | automatic | crystal | gamma } tpiba : read k-points in cartesian coordinates, in units of 2 pi/a (default) automatic: automatically generated uniform grid of k-points next card: nk1, nk2, nk3, k1, k2, k3 generates ( nk1, nk2, nk3 ) grid with ( k1, k2, k3 ) offset nk1, nk2, nk3 as in Monkhorst-Pack grids k1, k2, k3 must be 0 ( no offset ) or 1 ( grid displaced by half a grid step in the corresponding direction ) BEWARE: only grids having the full symmetry of the crystal work with tetrahedra. Some grids with offset may not work. crystal : read k-points in crystal coordinates, i.e. in relative coordinates of the reciprocal lattive vectors gamma : use k = 0 (no need to list k-point specifications after card) In this case wavefunctions can be chosen as real, and specialized subroutines optimized for calculations at the gamma point are used (memory and cpu requirements are reduced by approximately one half) next card: nks number of supplied special points xk_x, xk_y, xk_z, wk special points in the irreducible Brillouin Zone of the lattice (with all symmetries) and weights (see the literature for lists of special points and the corresponding weights) If the symmetry is lower than the full symmetry of the lattice, additional points with appropriate weights are generated. In a non-scf calculation, weights do not affect the results. If you just need eigenvalues and eigenvectors (for instance, for a band-structure plot), weights can be set to any value (for instance all equal to 1) ------------------------------------------------------------------------------- CELL_PARAMETERS optional card, needed only if ibrav = 0 is specified, ignored otherwise Syntax: CELL_PARAMETERS { cubic | hexagonal } a(1,1) a(2,1) a(3,1) a(1,2) a(2,2) a(3,2) a(1,3) a(2,3) a(3,3) a(:,1) = crystal axis 1 alat units if celldm(1) was specified 2 2 a.u. if celldm(1)=0 3 3 Keyword "cubic" or "hexagonal" specify if you want to look for symmetries derived from the cubic symmetry group (default) or from the hexagonal symmetry group (assuming c axis as the z axis, a axis along the x axis) ------------------------------------------------------------------------------- CLIMBING_IMAGES optional card, needed only if calculation = 'neb' and CI_scheme = 'manual' Syntax: CLIMBING_IMAGES index1, index2, ..., indexN where index1, index2, ..., indexN are the indices of the images to which apply the Climbing Image procedure. If more than an image is specified they must be separated by a comma ------------------------------------------------------------------------------- CONSTRAINTS Ionic Constraints Syntax: CONSTRAINTS nconstr { constr_tol } constr_type(.) constr(1,.) constr(2,.) ... { constr_target(.) } Where: nconstr INTEGER, number of constraints constr_tol REAL, tolerance for keeping the constraints satisfied constr_type(.) CHARACTER, type of constrain : 'type_coord' : constraint on global coordination-number, i.e. the average number of atoms of type B surrounding the atoms of type A. The coordination is defined by using a Fermi-Dirac. (four indexes must be specified). 'atom_coord' : constraint on local coordination-number, i.e. the average number of atoms of type A surrounding a specific atom. The coordination is defined by using a Fermi-Dirac. (four indexes must be specified). 'distance' : constraint on interatomic distance (two atom indexes must be specified ). 'planar_angle' : constraint on planar angle (three atom indexes must be specified). 'torsional_angle' : constraint on torsional angle (four atom indexes must be specified). 'bennett_proj' : constraint on the projection onto a given direction of the vector defined by the position of one atom minus the center of mass of the others. ( Ch.H. Bennett in Diffusion in Solids, Recent Developments, Ed. by A.S. Nowick and J.J. Burton, New York 1975 ). constr(1,.) constr(2,.) ... REAL, these variables have different meanings for different constraint types: 'type_coord' : constr(1) is the first index of the atomic type involved constr(2) is the second index of the atomic type involved constr(3) is the cut-off radius for estimating the coordination constr(4) is a smoothing parameter 'atom_coord' : constr(1) is the atom index of the atom with constrained coordination constr(2) is the index of the atomic type involved in the coordination constr(3) is the cut-off radius for estimating the coordination constr(4) is a smoothing parameter 'distance' : atoms indices object of the constraint, as they appear in the 'ATOMIC_POSITION' CARD 'planar_angle', 'torsional_angle' : atoms indices object of the constraint, as they appear in the 'ATOMIC_POSITION' CARD (beware the order) 'bennett_proj' : constr(1) is the index of the atom whose position is constrained. constr(2:4) are the three coordinates of the vector that specifies the constraint direction. constr_target REAL, target for the constrain ( angles are specified in degrees ). This variable is optional. ------------------------------------------------------------------------------- COLLECTIVE_VARS Collective variables used for meta-dynamics calculations Syntax: COLLECTIVE_VARS ncolvars { tolerance } colvar_type(.) colvar(1,.) colvar(2,.) ... Where: ncolvars INTEGER, number of collective variables tolerance REAL, tolerance used for SHAKE. colvar_type(.) CHARACTER, type of collective variable : ... see the definition of constr_type in the CONSTRAINTS card. colvar(1,.) colvar(2,.) ... REAL, these variables have different meanings for different collective variable types. See the definition of constr in the CONSTRAINTS card. ------------------------------------------------------------------------------- ibrav is the structure index: ibrav structure celldm(2)-celldm(6) 0 "free", see above not used 1 cubic P (sc) not used 2 cubic F (fcc) not used 3 cubic I (bcc) not used 4 Hexagonal and Trigonal P celldm(3)=c/a 5 Trigonal R celldm(4)=cos(alpha) 6 Tetragonal P (st) celldm(3)=c/a 7 Tetragonal I (bct) celldm(3)=c/a 8 Orthorhombic P celldm(2)=b/a,celldm(3)=c/a 9 Orthorhombic base-centered(bco) celldm(2)=b/a,celldm(3)=c/a 10 Orthorhombic face-centered celldm(2)=b/a,celldm(3)=c/a 11 Orthorhombic body-centered celldm(2)=b/a,celldm(3)=c/a 12 Monoclinic P celldm(2)=b/a,celldm(3)=c/a, celldm(4)=cos(ab) 13 Monoclinic base-centered celldm(2)=b/a,celldm(3)=c/a, celldm(4)=cos(ab) 14 Triclinic celldm(2)= b/a, celldm(3)= c/a, celldm(4)= cos(bc), celldm(5)= cos(ac), celldm(6)= cos(ab) For P lattices: the special axis (c) is the z-axis, one basal-plane vector (a) is along x, the other basal-plane vector (b) is at angle gamma for monoclinic, at 120 degrees for trigonal and hexagonal lattices, at 90 degrees for cubic, tetragonal, orthorhombic lattices sc simple cubic ==================== a1 = a(1,0,0), a2 = a(0,1,0), a3 = a(0,0,1) fcc face centered cubic ==================== a1 = (a/2)(-1,0,1), a2 = (a/2)(0,1,1), a3 = (a/2)(-1,1,0). bcc body entered cubic ==================== a1 = (a/2)(1,1,1), a2 = (a/2)(-1,1,1), a3 = (a/2)(-1,-1,1). simple hexagonal and trigonal(p) ==================== a1 = a(1,0,0), a2 = a(-1/2,sqrt(3)/2,0), a3 = a(0,0,c/a). trigonal(r) =================== for these groups, the z-axis is chosen as the 3-fold axis, but the crystallographic vectors form a three-fold star around the z-axis, and the primitive cell is a simple rhombohedron. The crystallographic vectors are: a1 = a(tx,-ty,tz), a2 = a(0,2ty,tz), a3 = a(-tx,-ty,tz). where c=cos(alpha) is the cosine of the angle alpha between any pair of crystallographic vectors, tc, ty, tz are defined as tx=sqrt((1-c)/2), ty=sqrt((1-c)/6), tz=sqrt((1+2c)/3) simple tetragonal (p) ==================== a1 = a(1,0,0), a2 = a(0,1,0), a3 = a(0,0,c/a) body centered tetragonal (i) ================================ a1 = (a/2)(1,-1,c/a), a2 = (a/2)(1,1,c/a), a3 = (a/2)(-1,-1,c/a). simple orthorhombic (p) ============================= a1 = (a,0,0), a2 = (0,b,0), a3 = (0,0,c) bco base centered orthorhombic ============================= a1 = (a/2,b/2,0), a2 = (-a/2,b/2,0), a3 = (0,0,c) face centered orthorhombic ============================= a1 = (a/2,0,c/2), a2 = (a/2,b/2,0), a3 = (0,b/2,c/2) body centered orthorhombic ============================= a1 = (a/2,b/2,c/2), a2 = (-a/2,b/2,c/2), a3 = (-a/2,-b/2,c/2) monoclinic (p) ============================= a1 = (a,0,0), a2= (b*sin(gamma), b*cos(gamma), 0), a3 = (0, 0, c) where gamma is the angle between axis a and b base centered monoclinic ============================= a1 = ( a/2, 0, -c/2), a2 = (b*cos(gamma), b*sin(gamma), 0), a3 = ( a/2, 0, c/2), where gamma is the angle between axis a and b triclinic ============================= a1 = (a, 0, 0), a2 = (b*cos(gamma), b*sin(gamma), 0) a3 = (c*cos(beta), c*(cos(alpha)-cos(beta)cos(gamma))/sin(gamma), c*sqrt( 1 + 2*cos(alpha)cos(beta)cos(gamma) - cos(alpha)^2-cos(beta)^2-cos(gamma)^2 )/sin(gamma) ) where alpha is the angle between axis b and c beta is the angle between axis a and c gamma is the angle between axis a and b ----------------------------------------------------------------------