Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1# Copyright (C) 2008 CSC - Scientific Computing Ltd. 

2"""This module defines an ASE interface to VASP. 

3 

4Developed on the basis of modules by Jussi Enkovaara and John 

5Kitchin. The path of the directory containing the pseudopotential 

6directories (potpaw,potpaw_GGA, potpaw_PBE, ...) should be set 

7by the environmental flag $VASP_PP_PATH. 

8 

9The user should also set the environmental flag $VASP_SCRIPT pointing 

10to a python script looking something like:: 

11 

12 import os 

13 exitcode = os.system('vasp') 

14 

15Alternatively, user can set the environmental flag $VASP_COMMAND pointing 

16to the command use the launch vasp e.g. 'vasp' or 'mpirun -n 16 vasp' 

17 

18http://cms.mpi.univie.ac.at/vasp/ 

19""" 

20 

21import os 

22import warnings 

23import shutil 

24from os.path import join, isfile, islink 

25from typing import List, Sequence, Tuple 

26 

27import numpy as np 

28 

29import ase 

30from ase.calculators.calculator import kpts2ndarray 

31from ase.calculators.vasp.setups import get_default_setups 

32 

33 

34def format_kpoints(kpts, atoms, reciprocal=False, gamma=False): 

35 tokens = [] 

36 append = tokens.append 

37 

38 append('KPOINTS created by Atomic Simulation Environment\n') 

39 

40 if isinstance(kpts, dict): 

41 kpts = kpts2ndarray(kpts, atoms=atoms) 

42 reciprocal = True 

43 

44 shape = np.array(kpts).shape 

45 

46 # Wrap scalar in list if necessary 

47 if shape == (): 

48 kpts = [kpts] 

49 shape = (1, ) 

50 

51 if len(shape) == 1: 

52 append('0\n') 

53 if shape == (1, ): 

54 append('Auto\n') 

55 elif gamma: 

56 append('Gamma\n') 

57 else: 

58 append('Monkhorst-Pack\n') 

59 append(' '.join(f'{kpt:d}' for kpt in kpts)) 

60 append('\n0 0 0\n') 

61 elif len(shape) == 2: 

62 append('%i \n' % (len(kpts))) 

63 if reciprocal: 

64 append('Reciprocal\n') 

65 else: 

66 append('Cartesian\n') 

67 for n in range(len(kpts)): 

68 [append('%f ' % kpt) for kpt in kpts[n]] 

69 if shape[1] == 4: 

70 append('\n') 

71 elif shape[1] == 3: 

72 append('1.0 \n') 

73 return ''.join(tokens) 

74 

75 

76# Parameters that can be set in INCAR. The values which are None 

77# are not written and default parameters of VASP are used for them. 

78 

79float_keys = [ 

80 'aexx', # Fraction of exact/DFT exchange 

81 'aggac', # Fraction of gradient correction to correlation 

82 'aggax', # Fraction of gradient correction to exchange 

83 'aldac', # Fraction of LDA correlation energy 

84 'amin', # 

85 'amix', # 

86 'amix_mag', # 

87 'bmix', # tags for mixing 

88 'bmix_mag', # 

89 'cshift', # Complex shift for dielectric tensor calculation (LOPTICS) 

90 'deper', # relative stopping criterion for optimization of eigenvalue 

91 'ebreak', # absolute stopping criterion for optimization of eigenvalues 

92 # (EDIFF/N-BANDS/4) 

93 'efield', # applied electrostatic field 

94 'emax', # energy-range for DOSCAR file 

95 'emin', # 

96 'enaug', # Density cutoff 

97 'encut', # Planewave cutoff 

98 'encutgw', # energy cutoff for response function 

99 'encutfock', # FFT grid in the HF related routines 

100 'hfscreen', # attribute to change from PBE0 to HSE 

101 'kspacing', # determines the number of k-points if the KPOINTS 

102 # file is not present. KSPACING is the smallest 

103 # allowed spacing between k-points in units of 

104 # $\AA$^{-1}$. 

105 'potim', # time-step for ion-motion (fs) 

106 'nelect', # total number of electrons 

107 'param1', # Exchange parameter 

108 'param2', # Exchange parameter 

109 'pomass', # mass of ions in am 

110 'pstress', # add this stress to the stress tensor, and energy E = V * 

111 # pstress 

112 'sigma', # broadening in eV 

113 'smass', # Nose mass-parameter (am) 

114 'spring', # spring constant for NEB 

115 'time', # special control tag 

116 'weimin', # maximum weight for a band to be considered empty 

117 'zab_vdw', # vdW-DF parameter 

118 'zval', # ionic valence 

119 # The next keywords pertain to the VTST add-ons from Graeme Henkelman's 

120 # group at UT Austin 

121 'jacobian', # Weight of lattice to atomic motion 

122 'ddr', # (DdR) dimer separation 

123 'drotmax', # (DRotMax) number of rotation steps per translation step 

124 'dfnmin', # (DFNMin) rotational force below which dimer is not rotated 

125 'dfnmax', # (DFNMax) rotational force below which dimer rotation stops 

126 'sltol', # convergence ratio for minimum eigenvalue 

127 'sdr', # finite difference for setting up Lanczos matrix and step 

128 # size when translating 

129 'maxmove', # Max step for translation for IOPT > 0 

130 'invcurv', # Initial curvature for LBFGS (IOPT = 1) 

131 'timestep', # Dynamical timestep for IOPT = 3 and IOPT = 7 

132 'sdalpha', # Ratio between force and step size for IOPT = 4 

133 # The next keywords pertain to IOPT = 7 (i.e. FIRE) 

134 'ftimemax', # Max time step 

135 'ftimedec', # Factor to dec. dt 

136 'ftimeinc', # Factor to inc. dt 

137 'falpha', # Parameter for velocity damping 

138 'falphadec', # Factor to dec. alpha 

139 'clz', # electron count for core level shift 

140 'vdw_radius', # Cutoff radius for Grimme's DFT-D2 and DFT-D3 and 

141 # Tkatchenko and Scheffler's DFT-TS dispersion corrections 

142 'vdw_scaling', # Global scaling parameter for Grimme's DFT-D2 dispersion 

143 # correction 

144 'vdw_d', # Global damping parameter for Grimme's DFT-D2 and Tkatchenko 

145 # and Scheffler's DFT-TS dispersion corrections 

146 'vdw_cnradius', # Cutoff radius for calculating coordination number in 

147 # Grimme's DFT-D3 dispersion correction 

148 'vdw_s6', # Damping parameter for Grimme's DFT-D2 and DFT-D3 and 

149 # Tkatchenko and Scheffler's DFT-TS dispersion corrections 

150 'vdw_s8', # Damping parameter for Grimme's DFT-D3 dispersion correction 

151 'vdw_sr', # Scaling parameter for Grimme's DFT-D2 and DFT-D3 and 

152 # Tkatchenko and Scheffler's DFT-TS dispersion correction 

153 'vdw_a1', # Damping parameter for Grimme's DFT-D3 dispersion correction 

154 'vdw_a2', # Damping parameter for Grimme's DFT-D3 dispersion correction 

155 'eb_k', # solvent permitivity in Vaspsol 

156 'tau', # surface tension parameter in Vaspsol 

157 'langevin_gamma_l', # Friction for lattice degrees of freedom 

158 'pmass', # Mass for latice degrees of freedom 

159 'bparam', # B parameter for nonlocal VV10 vdW functional 

160 'cparam', # C parameter for nonlocal VV10 vdW functional 

161 'aldax', # Fraction of LDA exchange (for hybrid calculations) 

162 'tebeg', # 

163 'teend', # temperature during run 

164 'andersen_prob', # Probability of collision in Andersen thermostat 

165 'apaco', # Distance cutoff for pair correlation function calc. 

166 'auger_ecblo', # Undocumented parameter for Auger calculations 

167 'auger_edens', # Density of electrons in conduction band 

168 'auger_hdens', # Density of holes in valence band 

169 'auger_efermi', # Fixed Fermi level for Auger calculations 

170 'auger_evbhi', # Upper bound for valence band maximum 

171 'auger_ewidth', # Half-width of energy window function 

172 'auger_occ_fac_eeh', # Undocumented parameter for Auger calculations 

173 'auger_occ_fac_ehh', # Undocumented parameter for Auger calculations 

174 'auger_temp', # Temperature for Auger calculation 

175 'dq', # Finite difference displacement magnitude (NMR) 

176 'avgap', # Average gap (Model GW) 

177 'ch_sigma', # Broadening of the core electron absorption spectrum 

178 'bpotim', # Undocumented Bond-Boost parameter (GH patches) 

179 'qrr', # Undocumented Bond-Boost parameter (GH patches) 

180 'prr', # Undocumented Bond-Boost parameter (GH patches) 

181 'rcut', # Undocumented Bond-Boost parameter (GH patches) 

182 'dvmax', # Undocumented Bond-Boost parameter (GH patches) 

183 'bfgsinvcurv', # Initial curvature for BFGS (GH patches) 

184 'damping', # Damping parameter for LBFGS (GH patches) 

185 'efirst', # Energy of first NEB image (GH patches) 

186 'elast', # Energy of final NEB image (GH patches) 

187 'fmagval', # Force magnitude convergence criterion (GH patches) 

188 'cmbj', # Modified Becke-Johnson MetaGGA c-parameter 

189 'cmbja', # Modified Becke-Johnson MetaGGA alpha-parameter 

190 'cmbjb', # Modified Becke-Johnson MetaGGA beta-parameter 

191 'sigma_nc_k', # Width of ion gaussians (VASPsol) 

192 'sigma_k', # Width of dielectric cavidty (VASPsol) 

193 'nc_k', # Cavity turn-on density (VASPsol) 

194 'lambda_d_k', # Debye screening length (VASPsol) 

195 'ediffsol', # Tolerance for solvation convergence (VASPsol) 

196 'deg_threshold', # Degeneracy threshold 

197 'omegamin', # Minimum frequency for dense freq. grid 

198 'omegamax', # Maximum frequency for dense freq. grid 

199 'rtime', # Undocumented parameter 

200 'wplasma', # Undocumented parameter 

201 'wplasmai', # Undocumented parameter 

202 'dfield', # Undocumented parameter 

203 'omegatl', # Maximum frequency for coarse freq. grid 

204 'encutgwsoft', # Soft energy cutoff for response kernel 

205 'encutlf', # Undocumented parameter 

206 'scissor', # Scissor correction for GW/BSE calcs 

207 'dimer_dist', # Distance between dimer images 

208 'step_size', # Step size for finite difference in dimer calculation 

209 'step_max', # Maximum step size for dimer calculation 

210 'minrot', # Minimum rotation allowed in dimer calculation 

211 'dummy_mass', # Mass of dummy atom(s?) 

212 'shaketol', # Tolerance for SHAKE algorithm 

213 'shaketolsoft', # Soft tolerance for SHAKE algorithm 

214 'shakesca', # Scaling of each step taken in SHAKE algorithm 

215 'hills_stride', # Undocumented metadynamics parameter 

216 'hills_h', # Height (in eV) of gaussian bias for metadynamics 

217 'hills_w', # Width of gaussian bias for metadynamics 

218 'hills_k', # Force constant coupling dummy&real for metadynamics 

219 'hills_m', # Mass of dummy particle for use in metadynamics 

220 'hills_temperature', # Temp. of dummy particle for metadynamics 

221 'hills_andersen_prob', # Probability of thermostat coll. for metadynamics 

222 'hills_sqq', # Nose-hoover particle mass for metadynamics 

223 'dvvdelta0', # Undocumented parameter 

224 'dvvvnorm0', # Undocumented parameter 

225 'dvvminpotim', # Undocumented parameter 

226 'dvvmaxpotim', # Undocumented parameter 

227 'efermi', # Undocumented parameter 

228 'enchg', # Undocumented charge fitting parameter 

229 'tau0', # Undocumented charge fitting parameter 

230 'encut4o', # Cutoff energy for 4-center integrals (HF) 

231 'param3', # Undocumented HF parameter 

232 'model_eps0', # Undocumented HF parameter 

233 'model_alpha', # Undocumented HF parameter 

234 'qmaxfockae', # Undocumented HF parameter 

235 'hfscreenc', # Range-separated screening length for correlations 

236 'hfrcut', # Cutoff radius for HF potential kernel 

237 'encutae', # Undocumented parameter for all-electron density calc. 

238 'encutsubrotscf', # Undocumented subspace rotation SCF parameter 

239 'enini', # Cutoff energy for wavefunctions (?) 

240 'wc', # Undocumented mixing parameter 

241 'enmax', # Cutoff energy for wavefunctions (?) 

242 'scalee', # Undocumented parameter 

243 'eref', # Reference energy 

244 'epsilon', # Dielectric constant of bulk charged cells 

245 'rcmix', # Mixing parameter for core density in rel. core calcs. 

246 'esemicore', # Energetic lower bound for states considered "semicore" 

247 'external_pressure', # Pressure for NPT calcs., equivalent to PSTRESS 

248 'lj_radius', # Undocumented classical vdW parameter 

249 'lj_epsilon', # Undocumented classical vdW parameter 

250 'lj_sigma', # Undocumented classical vdW parameter 

251 'mbd_beta', # TS MBD vdW correction damping parameter 

252 'scsrad', # Cutoff radius for dipole-dipole interaction tensor in SCS 

253 'hitoler', # Iterative Hirschfeld partitioning tolerance 

254 'lambda', # "Spring constant" for magmom constraint calcs. 

255 'kproj_threshold', # Threshold for k-point projection scheme 

256 'maxpwamp', # Undocumented HF parameter 

257 'vcutoff', # Undocumented parameter 

258 'mdtemp', # Temperature for AIMD 

259 'mdgamma', # Undocumented AIMD parameter 

260 'mdalpha', # Undocumented AIMD parameter 

261 'ofield_kappa', # Bias potential strength for interface pinning method 

262 'ofield_q6_near', # Steinhardt-Nelson Q6 parameters for interface pinning 

263 'ofield_q6_far', # Steinhardt-Nelson Q6 parameters for interface pinning 

264 'ofield_a', # Target order parameter for interface pinning method 

265 'pthreshold', # Don't print timings for routines faster than this value 

266 'qltol', # Eigenvalue tolerance for Lanczos iteration (instanton) 

267 'qdr', # Step size for building Lanczos matrix & CG (instanton) 

268 'qmaxmove', # Max step size (instanton) 

269 'qdt', # Timestep for quickmin minimization (instanton) 

270 'qtpz', # Temperature (instanton) 

271 'qftol', # Tolerance (instanton) 

272 'nupdown', # fix spin moment to specified value 

273] 

274 

275exp_keys = [ 

276 'ediff', # stopping-criterion for electronic upd. 

277 'ediffg', # stopping-criterion for ionic upd. 

278 'symprec', # precession in symmetry routines 

279 # The next keywords pertain to the VTST add-ons from Graeme Henkelman's 

280 # group at UT Austin 

281 'fdstep', # Finite diference step for IOPT = 1 or 2 

282] 

283 

284string_keys = [ 

285 'algo', # algorithm: Normal (Davidson) | Fast | Very_Fast (RMM-DIIS) 

286 'gga', # xc-type: PW PB LM or 91 (LDA if not set) 

287 'metagga', # 

288 'prec', # Precission of calculation (Low, Normal, Accurate) 

289 'system', # name of System 

290 'precfock', # FFT grid in the HF related routines 

291 'radeq', # Which type of radial equations to use for rel. core calcs. 

292 'localized_basis', # Basis to use in CRPA 

293 'proutine', # Select profiling routine 

294] 

295 

296int_keys = [ 

297 'ialgo', # algorithm: use only 8 (CG) or 48 (RMM-DIIS) 

298 'ibrion', # ionic relaxation: 0-MD 1-quasi-New 2-CG 

299 'icharg', # charge: 0-WAVECAR 1-CHGCAR 2-atom 10-const 

300 'idipol', # monopol/dipol and quadropole corrections 

301 'images', # number of images for NEB calculation 

302 'iniwav', # initial electr wf. : 0-lowe 1-rand 

303 'isif', # calculate stress and what to relax 

304 'ismear', # part. occupancies: -5 Blochl -4-tet -1-fermi 0-gaus >0 MP 

305 'ispin', # spin-polarized calculation 

306 'istart', # startjob: 0-new 1-cont 2-samecut 

307 'isym', # symmetry: 0-nonsym 1-usesym 2-usePAWsym 

308 'iwavpr', # prediction of wf.: 0-non 1-charg 2-wave 3-comb 

309 'kpar', # k-point parallelization paramater 

310 'ldauprint', # 0-silent, 1-occ. matrix written to OUTCAR, 2-1+pot. matrix 

311 # written 

312 'ldautype', # L(S)DA+U: 1-Liechtenstein 2-Dudarev 4-Liechtenstein(LDAU) 

313 'lmaxmix', # 

314 'lorbit', # create PROOUT 

315 'maxmix', # 

316 'ngx', # FFT mesh for wavefunctions, x 

317 'ngxf', # FFT mesh for charges x 

318 'ngy', # FFT mesh for wavefunctions, y 

319 'ngyf', # FFT mesh for charges y 

320 'ngz', # FFT mesh for wavefunctions, z 

321 'ngzf', # FFT mesh for charges z 

322 'nbands', # Number of bands 

323 'nblk', # blocking for some BLAS calls (Sec. 6.5) 

324 'nbmod', # specifies mode for partial charge calculation 

325 'nelm', # nr. of electronic steps (default 60) 

326 'nelmdl', # nr. of initial electronic steps 

327 'nelmgw', # nr. of self-consistency cycles for GW 

328 'nelmin', 

329 'nfree', # number of steps per DOF when calculting Hessian using 

330 # finite differences 

331 'nkred', # define sub grid of q-points for HF with 

332 # nkredx=nkredy=nkredz 

333 'nkredx', # define sub grid of q-points in x direction for HF 

334 'nkredy', # define sub grid of q-points in y direction for HF 

335 'nkredz', # define sub grid of q-points in z direction for HF 

336 'nomega', # number of frequency points 

337 'nomegar', # number of frequency points on real axis 

338 'npar', # parallelization over bands 

339 'nsim', # evaluate NSIM bands simultaneously if using RMM-DIIS 

340 'nsw', # number of steps for ionic upd. 

341 'nwrite', # verbosity write-flag (how much is written) 

342 'vdwgr', # extra keyword for Andris program 

343 'vdwrn', # extra keyword for Andris program 

344 'voskown', # use Vosko, Wilk, Nusair interpolation 

345 # The next keywords pertain to the VTST add-ons from Graeme Henkelman's 

346 # group at UT Austin 

347 'ichain', # Flag for controlling which method is being used (0=NEB, 

348 # 1=DynMat, 2=Dimer, 3=Lanczos) if ichain > 3, then both 

349 # IBRION and POTIM are automatically set in the INCAR file 

350 'iopt', # Controls which optimizer to use. for iopt > 0, ibrion = 3 

351 # and potim = 0.0 

352 'snl', # Maximum dimentionality of the Lanczos matrix 

353 'lbfgsmem', # Steps saved for inverse Hessian for IOPT = 1 (LBFGS) 

354 'fnmin', # Max iter. before adjusting dt and alpha for IOPT = 7 (FIRE) 

355 'icorelevel', # core level shifts 

356 'clnt', # species index 

357 'cln', # main quantum number of excited core electron 

358 'cll', # l quantum number of excited core electron 

359 'ivdw', # Choose which dispersion correction method to use 

360 'nbandsgw', # Number of bands for GW 

361 'nbandso', # Number of occupied bands for electron-hole treatment 

362 'nbandsv', # Number of virtual bands for electron-hole treatment 

363 'ncore', # Number of cores per band, equal to number of cores divided 

364 # by npar 

365 'mdalgo', # Determines which MD method of Tomas Bucko to use 

366 'nedos', # Number of grid points in DOS 

367 'turbo', # Ewald, 0 = Normal, 1 = PME 

368 'omegapar', # Number of groups for response function calc. 

369 # (Possibly Depricated) Number of groups in real time for 

370 # response function calc. 

371 'taupar', 

372 'ntaupar', # Number of groups in real time for response function calc. 

373 'antires', # How to treat antiresonant part of response function 

374 'magatom', # Index of atom at which to place magnetic field (NMR) 

375 'jatom', # Index of atom at which magnetic moment is evaluated (NMR) 

376 'ichibare', # chi_bare stencil size (NMR) 

377 'nbas', # Undocumented Bond-Boost parameter (GH patches) 

378 'rmds', # Undocumented Bond-Boost parameter (GH patches) 

379 'ilbfgsmem', # Number of histories to store for LBFGS (GH patches) 

380 'vcaimages', # Undocumented parameter (GH patches) 

381 'ntemper', # Undocumented subspace diagonalization param. (GH patches) 

382 'ncshmem', # Share memory between this many cores on each process 

383 'lmaxtau', # Undocumented MetaGGA parameter (prob. max ang.mom. for tau) 

384 'kinter', # Additional finer grid (?) 

385 'ibse', # Type of BSE calculation 

386 'nbseeig', # Number of BSE wfns to write 

387 'naturalo', # Use NATURALO (?) 

388 'nbandsexact', # Undocumented parameter 

389 'nbandsgwlow', # Number of bands for which shifts are calculated 

390 'nbandslf', # Number of bands included in local field effect calc. 

391 'omegagrid', # Undocumented parameter 

392 'telescope', # Undocumented parameter 

393 'maxmem', # Amount of memory to allocate per core in MB 

394 'nelmhf', # Number of iterations for HF part (GW) 

395 'dim', # Undocumented parameter 

396 'nkredlf', # Reduce k-points for local field effects 

397 'nkredlfx', # Reduce k-points for local field effects in X 

398 'nkredlfy', # Reduce k-points for local field effects in Y 

399 'nkredlfz', # Reduce k-points for local field effects in Z 

400 'lmaxmp2', # Undocumented parameter 

401 'switch', # Undocumented dimer parameter 

402 'findiff', # Use forward (1) or central (2) finite difference for dimer 

403 'engine', # Undocumented dimer parameter 

404 'restartcg', # Undocumented dimer parameter 

405 'thermostat', # Deprecated parameter for selecting MD method (use MDALGO) 

406 'scaling', # After how many steps velocities should be rescaled 

407 'shakemaxiter', # Maximum # of iterations in SHAKE algorithm 

408 'equi_regime', # Number of steps to equilibrate for 

409 'hills_bin', # Update metadynamics bias after this many steps 

410 'hills_maxstride', # Undocumented metadynamics parameter 

411 'dvvehistory', # Undocumented parameter 

412 'ipead', # Undocumented parameter 

413 'ngaus', # Undocumented charge fitting parameter 

414 'exxoep', # Undocumented HF parameter 

415 'fourorbit', # Undocumented HF parameter 

416 'model_gw', # Undocumented HF parameter 

417 'hflmax', # Maximum L quantum number for HF calculation 

418 'lmaxfock', # Maximum L quantum number for HF calc. (same as above) 

419 'lmaxfockae', # Undocumented HF parameter 

420 'nmaxfockae', # Undocumented HF parameter 

421 'nblock_fock', # Undocumented HF parameter 

422 'idiot', # Determines which warnings/errors to print 

423 'nrmm', # Number of RMM-DIIS iterations 

424 'mremove', # Undocumented mixing parameter 

425 'inimix', # Undocumented mixing parameter 

426 'mixpre', # Undocumented mixing parameter 

427 'nelmall', # Undocumented parameter 

428 'nblock', # How frequently to write data 

429 'kblock', # How frequently to write data 

430 'npaco', # Undocumented pair correlation function parameter 

431 'lmaxpaw', # Max L quantum number for on-site charge expansion 

432 'irestart', # Undocumented parameter 

433 'nreboot', # Undocumented parameter 

434 'nmin', # Undocumented parameter 

435 'nlspline', # Undocumented parameter 

436 'ispecial', # "Select undocumented and unsupported special features" 

437 'rcrep', # Number of steps between printing relaxed core info 

438 'rcndl', # Wait this many steps before updating core density 

439 'rcstrd', # Relax core density after this many SCF steps 

440 'vdw_idampf', # Select type of damping function for TS vdW 

441 'i_constrained_m', # Select type of magmom. constraint to use 

442 'igpar', # "G parallel" direction for Berry phase calculation 

443 'nppstr', # Number of kpts in "igpar' direction for Berry phase calc. 

444 'nbands_out', # Undocumented QP parameter 

445 'kpts_out', # Undocumented QP parameter 

446 'isp_out', # Undocumented QP parameter 

447 'nomega_out', # Undocumented QP parameter 

448 'maxiter_ft', # Max iterations for sloppy Remez algorithm 

449 'nmaxalt', # Max sample points for alternant in Remez algorithms 

450 'itmaxlsq', # Max iterations in LSQ search algorithm 

451 'ndatalsq', # Number of sample points for LSQ search algorithm 

452 'ncore_in_image1', # Undocumented parameter 

453 'kimages', # Undocumented parameter 

454 'ncores_per_band', # Undocumented parameter 

455 'maxlie', # Max iterations in CRPA diagonalization routine 

456 'ncrpalow', # Undocumented CRPA parameter 

457 'ncrpahigh', # Undocumented CRPA parameter 

458 'nwlow', # Undocumented parameter 

459 'nwhigh', # Undocumented parameter 

460 'nkopt', # Number of k-points to include in Optics calculation 

461 'nkoffopt', # K-point "counter offset" for Optics 

462 'nbvalopt', # Number of valence bands to write in OPTICS file 

463 'nbconopt', # Number of conduction bands to write in OPTICS file 

464 'ch_nedos', # Number dielectric function calculation grid points for XAS 

465 'plevel', # No timings for routines with "level" higher than this 

466 'qnl', # Lanczos matrix size (instanton) 

467] 

468 

469bool_keys = [ 

470 'addgrid', # finer grid for augmentation charge density 

471 'kgamma', # The generated kpoint grid (from KSPACING) is either 

472 # centred at the $\Gamma$ 

473 # point (e.g. includes the $\Gamma$ point) 

474 # (KGAMMA=.TRUE.) 

475 'laechg', # write AECCAR0/AECCAR1/AECCAR2 

476 'lasph', # non-spherical contributions to XC energy (and pot for 

477 # VASP.5.X) 

478 'lasync', # overlap communcation with calculations 

479 'lcharg', # 

480 'lcorr', # Harris-correction to forces 

481 'ldau', # L(S)DA+U 

482 'ldiag', # algorithm: perform sub space rotation 

483 'ldipol', # potential correction mode 

484 'lelf', # create ELFCAR 

485 'lepsilon', # enables to calculate and to print the BEC tensors 

486 'lhfcalc', # switch to turn on Hartree Fock calculations 

487 'loptics', # calculate the frequency dependent dielectric matrix 

488 'lpard', # evaluate partial (band and/or k-point) decomposed charge 

489 # density 

490 'lplane', # parallelisation over the FFT grid 

491 'lscalapack', # switch off scaLAPACK 

492 'lscalu', # switch of LU decomposition 

493 'lsepb', # write out partial charge of each band separately? 

494 'lsepk', # write out partial charge of each k-point separately? 

495 'lthomas', # 

496 'luse_vdw', # Invoke vdW-DF implementation by Klimes et. al 

497 'lvdw', # Invoke DFT-D2 method of Grimme 

498 'lvhar', # write Hartree potential to LOCPOT (vasp 5.x) 

499 'lvtot', # create WAVECAR/CHGCAR/LOCPOT 

500 'lwave', # 

501 # The next keywords pertain to the VTST add-ons from Graeme Henkelman's 

502 # group at UT Austin 

503 'lclimb', # Turn on CI-NEB 

504 'ltangentold', # Old central difference tangent 

505 'ldneb', # Turn on modified double nudging 

506 'lnebcell', # Turn on SS-NEB 

507 'lglobal', # Optmize NEB globally for LBFGS (IOPT = 1) 

508 'llineopt', # Use force based line minimizer for translation (IOPT = 1) 

509 'lbeefens', # Switch on print of BEE energy contributions in OUTCAR 

510 'lbeefbas', # Switch off print of all BEEs in OUTCAR 

511 'lcalcpol', # macroscopic polarization (vasp5.2). 'lcalceps' 

512 'lcalceps', # Macroscopic dielectric properties and Born effective charge 

513 # tensors (vasp 5.2) 

514 'lvdw', # Turns on dispersion correction 

515 'lvdw_ewald', # Turns on Ewald summation for Grimme's DFT-D2 and 

516 # Tkatchenko and Scheffler's DFT-TS dispersion correction 

517 'lspectral', # Use the spectral method to calculate independent particle 

518 # polarizability 

519 'lrpa', # Include local field effects on the Hartree level only 

520 'lwannier90', # Switches on the interface between VASP and WANNIER90 

521 'lsorbit', # Enable spin-orbit coupling 

522 'lsol', # turn on solvation for Vaspsol 

523 'lautoscale', # automatically calculate inverse curvature for VTST LBFGS 

524 'interactive', # Enables interactive calculation for VaspInteractive 

525 'lauger', # Perform Auger calculation (Auger) 

526 'lauger_eeh', # Calculate EEH processes (Auger) 

527 'lauger_ehh', # Calculate EHH processes (Auger) 

528 'lauger_collect', # Collect wfns before looping over k-points (Auger) 

529 'lauger_dhdk', # Auto-determine E. window width from E. derivs. (Auger) 

530 'lauger_jit', # Distribute wavefunctions for k1-k4 (Auger) 

531 'orbitalmag', # Enable orbital magnetization (NMR) 

532 'lchimag', # Use linear response for shielding tensor (NMR) 

533 'lwrtcur', # Write response of current to mag. field to file (NMR) 

534 'lnmr_sym_red', # Reduce symmetry for finite difference (NMR) 

535 'lzora', # Use ZORA approximation in linear-response NMR (NMR) 

536 'lbone', # Use B-component in AE one-center terms for LR NMR (NMR) 

537 'lmagbloch', # Use Bloch summations to obtain orbital magnetization (NMR) 

538 'lgauge', # Use gauge transformation for zero moment terms (NMR) 

539 'lbfconst', # Use constant B-field with sawtooth vector potential (NMR) 

540 'nucind', # Use nuclear independent calculation (NMR) 

541 'lnicsall', # Use all grid points for 'nucind' calculation (NMR) 

542 'llraug', # Use two-center corrections for induced B-field (NMR) 

543 'lbbm', # Undocumented Bond-Boost parameter (GH patches) 

544 'lnoncollinear', # Do non-collinear spin polarized calculation 

545 'bfgsdfp', # Undocumented BFGS parameter (GH patches) 

546 'linemin', # Use line minimization (GH patches) 

547 'ldneborg', # Undocumented NEB parameter (GH patches) 

548 'dseed', # Undocumented dimer parameter (GH patches) 

549 'linteract', # Undocumented parameter (GH patches) 

550 'lmpmd', # Undocumented parameter (GH patches) 

551 'ltwodim', # Makes stress tensor two-dimensional (GH patches) 

552 'fmagflag', # Use force magnitude as convergence criterion (GH patches) 

553 'ltemper', # Use subspace diagonalization (?) (GH patches) 

554 'qmflag', # Undocumented FIRE parameter (GH patches) 

555 'lmixtau', # Undocumented MetaGGA parameter 

556 'ljdftx', # Undocumented VASPsol parameter (VASPsol) 

557 'lrhob', # Write the bound charge density (VASPsol) 

558 'lrhoion', # Write the ionic charge density (VASPsol) 

559 'lnabla', # Undocumented parameter 

560 'linterfast', # Interpolate in K using linear response routines 

561 'lvel', # Undocumented parameter 

562 'lrpaforce', # Calculate RPA forces 

563 'lhartree', # Use IP approx. in BSE (testing only) 

564 'ladder', # Use ladder diagrams 

565 'lfxc', # Use approximate ladder diagrams 

566 'lrsrpa', # Undocumented parameter 

567 'lsingles', # Calculate HF singles 

568 'lfermigw', # Iterate Fermi level 

569 'ltcte', # Undocumented parameter 

570 'ltete', # Undocumented parameter 

571 'ltriplet', # Undocumented parameter 

572 'lfxceps', # Undocumented parameter 

573 'lfxheg', # Undocumented parameter 

574 'l2order', # Undocumented parameter 

575 'lmp2lt', # Undocumented parameter 

576 'lgwlf', # Undocumented parameter 

577 'lusew', # Undocumented parameter 

578 'selfenergy', # Undocumented parameter 

579 'oddonlygw', # Avoid gamma point in response function calc. 

580 'evenonlygw', # Avoid even points in response function calc. 

581 'lspectralgw', # More accurate self-energy calculation 

582 'ch_lspec', # Calculate matrix elements btw. core and conduction states 

583 'fletcher_reeves', # Undocumented dimer parameter 

584 'lidm_selective', # Undocumented dimer parameter 

585 'lblueout', # Write output of blue-moon algorithm 

586 'hills_variable_w', # Enable variable-width metadynamics bias 

587 'dvvminus', # Undocumented parameter 

588 'lpead', # Calculate cell-periodic orbital derivs. using finite diff. 

589 'skip_edotp', # Skip updating elec. polarization during scf 

590 'skip_scf', # Skip calculation w/ local field effects 

591 'lchgfit', # Turn on charge fitting 

592 'lgausrc', # Undocumented charge fitting parameter 

593 'lstockholder', # Enable ISA charge fitting (?) 

594 'lsymgrad', # Restore symmetry of gradient (HF) 

595 'lhfone', # Calculate one-center terms (HF) 

596 'lrscor', # Include long-range correlation (HF) 

597 'lrhfcalc', # Include long-range HF (HF) 

598 'lmodelhf', # Model HF calculation (HF) 

599 'shiftred', # Undocumented HF parameter 

600 'hfkident', # Undocumented HF parameter 

601 'oddonly', # Undocumented HF parameter 

602 'evenonly', # Undocumented HF parameter 

603 'lfockaedft', # Undocumented HF parameter 

604 'lsubrot', # Enable subspace rotation diagonalization 

605 'mixfirst', # Mix before diagonalization 

606 'lvcader', # Calculate derivs. w.r.t. VCA parameters 

607 'lcompat', # Enable "full compatibility" 

608 'lmusic', # "Joke" parameter 

609 'ldownsample', # Downsample WAVECAR to fewer k-points 

610 'lscaaware', # Disable ScaLAPACK for some things but not all 

611 'lorbitalreal', # Undocumented parameter 

612 'lmetagga', # Undocumented parameter 

613 'lspiral', # Undocumented parameter 

614 'lzeroz', # Undocumented parameter 

615 'lmono', # Enable "monopole" corrections 

616 'lrelcore', # Perform relaxed core calculation 

617 'lmimicfc', # Mimic frozen-core calcs. for relaxed core calcs. 

618 'lmatchrw', # Match PS partial waves at RWIGS? (otherwise PAW cutoff) 

619 'ladaptelin', # Linearize core state energies to avoid divergences 

620 'lonlysemicore', # Only linearize semi-core state energies 

621 'gga_compat', # Enable backwards-compatible symmetrization of GGA derivs. 

622 'lrelvol', # Undocumented classical vdW parameter 

623 'lj_only', # Undocumented classical vdW parameter 

624 'lvdwscs', # Include self-consistent screening in TS vdW correction 

625 'lcfdm', # Use coupled fluctuating dipoles model for TS vdW 

626 'lvdw_sametype', # Include interactions between atoms of the same type 

627 'lrescaler0', # Rescale damping parameters in SCS vdW correction 

628 'lscsgrad', # Calculate gradients for TS+SCS vdW correction energies 

629 'lvdwexpansion', # Write 2-6 body contribs. to MBD vdW correction energy 

630 'lvdw_relvolone', # Undocumented classical vdW parameter 

631 'lberry', # Enable Berry-phase calculation 

632 'lpade_fit', # Undocumented QP parameter 

633 'lkproj', # Enable projection onto k-points 

634 'l_wr_moments', # Undocumented parameter 

635 'l_wr_density', # Undocumented parameter 

636 'lkotani', # Undocumented parameter 

637 'ldyson', # Undocumented parameter 

638 'laddherm', # Undocumented parameter 

639 'lcrpaplot', # Plot bands used in CRPA response func. calc. 

640 'lplotdis', # Plot disentangled bands in CRPA response func. calc. 

641 'ldisentangle', # Disentangle bands in CRPA 

642 'lweighted', # "Weighted" CRPA approach 

643 'luseorth_lcaos', # Use orthogonalized LCAOs in CRPA 

644 'lfrpa', # Use full RPA in CRPA 

645 'lregularize', # Regularize projectors in CRPA 

646 'ldrude', # Include Drude term in CRPA 

647 'ldmatrix', # Undocumented parameter 

648 'lefg', # Calculate electric field gradient at atomic nuclei 

649 'lhyperfine', # Enable Hyperfine calculation 

650 'lwannier', # Enable Wannier interface 

651 'localize', # Undocumented Wannier parameter 

652 'lintpol_wpot', # Interpolate WPOT for Wannier 

653 'lintpol_orb', # Interpolate orbitals for Wannier 

654 'lintpol_kpath', # Interpolate bandstructure on given kpath for Wannier 

655 'lintpol_kpath_orb', # Interpolate orbitals on given kpath for Wannier 

656 'lread_eigenvalues', # Use Eigenvalues from EIGENVALUES.INT file 

657 'lintpol_velocity', # Interpolate electron velocity for Wannier 

658 'lintpol_conductivity', # Interpolate conductivity for Wannier 

659 'lwannierinterpol', # Undocumented Wannier parameter 

660 'wanproj', # Undocumented Wannier parameter 

661 'lorbmom', # Undocumented LDA+U parameter 

662 'lwannier90_run', # Undocumented WANNIER90 parameter 

663 'lwrite_wanproj', # Write UWAN files for WANNIER90 

664 'lwrite_unk', # Write UNK files for WANNIER90 

665 'lwrite_mmn_amn', # Write MMN and AMN files for WANNIER90 

666 'lread_amn', # Read AMN files instead of recomputing (WANNIER90) 

667 'lrhfatm', # Undocumented HF parameter 

668 'lvpot', # Calculate unscreened potential 

669 'lwpot', # Calculate screened potential 

670 'lwswq', # Undocumented parameter 

671 'pflat', # Only print "flat" timings to OUTCAR 

672 'qifcg', # Use CG instead of quickmin (instanton) 

673 'qdo_ins', # Find instanton 

674 'qdo_pre', # Calculate prefactor (instanton) 

675 # The next keyword pertains to the periodic NBO code of JR Schmidt's group 

676 # at UW-Madison (https://github.com/jrschmidt2/periodic-NBO) 

677 'lnbo', # Enable NBO analysis 

678] 

679 

680list_int_keys = [ 

681 'iband', # bands to calculate partial charge for 

682 'kpuse', # k-point to calculate partial charge for 

683 'ldaul', # DFT+U parameters, overruled by dict key 'ldau_luj' 

684 'random_seed', # List of ints used to seed RNG for advanced MD routines 

685 # (Bucko) 

686 'auger_bmin_eeh', # 4 ints | Various undocumented parameters for Auger 

687 'auger_bmax_eeh', # 4 ints | calculations 

688 'auger_bmin_ehh', # 4 ints | 

689 'auger_bmax_ehh', # 4 ints | 

690 'balist', # nbas ints | Undocumented Bond-Boost parameter (GH patches) 

691 'kpoint_bse', # 4 ints | Undocumented parameter 

692 'nsubsys', # <=3 ints | Last atom # for each of up to 3 thermostats 

693 'vdw_refstate', # ntyp ints | Undocumented classical vdW parameter 

694 'vdw_mbd_size', # 3 ints | Supercell size for TS MBD vdW correction 

695 'nbands_index', # nbands_out ints | Undocumented QP parameter 

696 'kpts_index', # kpts_out ints | Undocumented QP parameter 

697 'isp_index', # isp_out ints | Undocumented QP parameter 

698 'nomega_index', # nomega_out ints | Undocumented QP parameter 

699 'ntarget_states', # nbands ints | Undocumented CRPA parameter 

700 'wanproj_i', # nions ints | Undocumented Wannier parameter 

701 'wanproj_l', # ? ints | Undocumented Wannier parameter 

702] 

703 

704list_bool_keys = [ 

705 'lattice_constraints', # 3 bools | Undocumented advanced MD parameter 

706 'lrctype', # ntyp bools | Enable relaxed-core calc. for these atoms 

707 'lvdw_onecell', # 3 bools | Enable periodicity in A, B, C vector for vdW 

708] 

709 

710list_float_keys = [ 

711 'dipol', # center of cell for dipol 

712 'eint', # energy range to calculate partial charge for 

713 'ferwe', # Fixed band occupation (spin-paired) 

714 'ferdo', # Fixed band occupation (spin-plarized) 

715 'magmom', # initial magnetic moments 

716 'ropt', # number of grid points for non-local proj in real space 

717 'rwigs', # Wigner-Seitz radii 

718 'ldauu', # ldau parameters, has potential to redundant w.r.t. dict 

719 'ldauj', # key 'ldau_luj', but 'ldau_luj' can't be read direct from 

720 # the INCAR (since it needs to know information about atomic 

721 # species. In case of conflict 'ldau_luj' gets written out 

722 # when a calculation is set up 

723 'vdw_c6', # List of floats of C6 parameters (J nm^6 mol^-1) for each 

724 # species (DFT-D2 and DFT-TS) 

725 'vdw_c6au', # List of floats of C6 parameters (a.u.) for each species 

726 # (DFT-TS) 

727 'vdw_r0', # List of floats of R0 parameters (angstroms) for each 

728 # species (DFT-D2 and DFT-TS) 

729 'vdw_r0au', # List of floats of R0 parameters (a.u.) for each species 

730 # (DFT-TS) 

731 'vdw_alpha', # List of floats of free-atomic polarizabilities for each 

732 # species (DFT-TS) 

733 'langevin_gamma', # List of floats for langevin friction coefficients 

734 'auger_emin_eeh', # 4 floats | Various undocumented parameters for Auger 

735 'auger_emax_eeh', # 4 floats | calculations 

736 'auger_emin_ehh', # 4 floats | 

737 'auger_emax_ehh', # 4 floats | 

738 'avecconst', # 3 floats | magnitude of magnetic moment (NMR) 

739 'magdipol', # 3 floats | magnitude of magnetic dipole (NMR) 

740 'bconst', # 3 floats | magnitude of constant magnetic field (NMR) 

741 'magpos', # 3 floats | position for magnetic moment w/ 'nucind' (NMR) 

742 'bext', # 3 floats | Undocumented (probably external magnetic field) 

743 'core_c', # ntyp floats | pseudo-core charge magnitude (VASPsol) 

744 'sigma_rc_k', # ntyp floats | width of pseudo-core gaussians (VASPsol) 

745 'darwinr', # ntypd (?) floats | Undocumented parameter 

746 'darwinv', # ntypd (?) floats | Undocumented parameter 

747 'dummy_k', # ? floats | Force const. connecting dummy atoms to sys. 

748 'dummy_r0', # ? floats | Minimum dist., ang., etc. for dummy atom DOFs 

749 'dummy_positions', # 3 floats | Position of dummy atom(s?) 

750 'psubsys', # <=3 floats | Coll. prob. for each of up to 3 thermostats 

751 'tsubsys', # <=3 floats | Temp. for each of up to 3 thermostats 

752 'increm', # ? floats | Undocumented advanced MD parameter 

753 'value_min', # ? floats | Undocumented advanced MD parameter 

754 'value_max', # ? floats | Undocumented advanced MD parameter 

755 'hills_position', # ? floats | Dummy particle(s) pos. for metadynamics 

756 'hills_velocity', # ? floats | Dummy particle(s) vel. for metadynamics 

757 'spring_k', # ? floats | Spring constant for harmonic constraints 

758 'spring_r0', # ? floats | Spring minima for harmonic constraints 

759 'spring_v0', # ? floats | Initial velocity of harmonic constraints 

760 'hills_wall_lower', # ? floats | Undocumented metadynamics parameter 

761 'hills_wall_upper', # ? floats | Undocumented metadynamics parameter 

762 'efield_pead', # 3 floats | homogeneous electric field for PEAD calc. 

763 'zct', # ? floats | Undocumented charge fitting parameter 

764 'rgaus', # ? floats | Undocumented charge fitting parameter 

765 'hfalpha', # 10 floats | Undocumented HF parameter 

766 'mcalpha', # 10 floats | Undocumented HF parameter 

767 'saxis', # 3 floats | Coordinate for collinear spin calculations 

768 'vca', # ? floats | Atom weight for VCA calculations 

769 'stm', # 7 floats | "range for STM data" 

770 'qspiral', # 3 floats | Undocumented parameter 

771 'external_stress', # 6 floats | Target stress (adds w/ external_pressure) 

772 'm_constr', # 3*nions floats | Local magmom assigned to each spin DOF 

773 'quad_efg', # ntyp floats | Nuclear quadrupole moments 

774 'ngyromag', # ntyp floats | Nuclear gyromagnetic ratios 

775 'rcrhocut', # ntyp floats | Core density cutoff rad. for HF relcore calc 

776 'ofield_k', # 3 floats | Undocumented parameter 

777 'paripot', # ? floats | Undocumented parameter 

778 'smearings', # ? floats | ismear,sigma smearing params to loop over 

779 'wanproj_e', # 2 floats | Undocumented Wannier parameter 

780] 

781 

782special_keys = [ 

783 'lreal', # non-local projectors in real space 

784] 

785 

786dict_keys = [ 

787 'ldau_luj', # dictionary with L(S)DA+U parameters, e.g. {'Fe':{'L':2, 

788 # 'U':4.0, 'J':0.9}, ...} 

789] 

790 

791keys: List[str] = [ 

792 # 'NBLOCK' and KBLOCK inner block; outer block 

793 # 'NPACO' and APACO distance and nr. of slots for P.C. 

794 # 'WEIMIN, EBREAK, DEPER special control tags 

795] 

796 

797 

798class GenerateVaspInput: 

799 # Parameters corresponding to 'xc' settings. This may be modified 

800 # by the user in-between loading calculators.vasp submodule and 

801 # instantiating the calculator object with calculators.vasp.Vasp() 

802 xc_defaults = { 

803 'lda': { 

804 'pp': 'LDA' 

805 }, 

806 # GGAs 

807 'blyp': { # https://www.vasp.at/forum/viewtopic.php?p=17234 

808 'pp': 'PBE', 

809 'gga': 'B5', 

810 'aldax': 1.00, 

811 'aggax': 1.00, 

812 'aggac': 1.00, 

813 'aldac': 0.00 

814 }, 

815 'pw91': { 

816 'pp': 'PW91', 

817 'gga': '91' 

818 }, 

819 'pbe': { 

820 'pp': 'PBE', 

821 'gga': 'PE' 

822 }, 

823 'pbesol': { 

824 'gga': 'PS' 

825 }, 

826 'revpbe': { 

827 'gga': 'RE' 

828 }, 

829 'rpbe': { 

830 'gga': 'RP' 

831 }, 

832 'am05': { 

833 'gga': 'AM' 

834 }, 

835 # Meta-GGAs 

836 'tpss': { 

837 'metagga': 'TPSS' 

838 }, 

839 'revtpss': { 

840 'metagga': 'RTPSS' 

841 }, 

842 'm06l': { 

843 'metagga': 'M06L' 

844 }, 

845 'ms0': { 

846 'metagga': 'MS0' 

847 }, 

848 'ms1': { 

849 'metagga': 'MS1' 

850 }, 

851 'ms2': { 

852 'metagga': 'MS2' 

853 }, 

854 'scan': { 

855 'metagga': 'SCAN' 

856 }, 

857 'rscan': { 

858 'metagga': 'RSCAN' 

859 }, 

860 'r2scan': { 

861 'metagga': 'R2SCAN' 

862 }, 

863 'scan-rvv10': { 

864 'metagga': 'SCAN', 

865 'luse_vdw': True, 

866 'bparam': 15.7 

867 }, 

868 'mbj': { 

869 # Modified Becke-Johnson 

870 'metagga': 'MBJ', 

871 }, 

872 'tb09': { 

873 # Alias for MBJ 

874 'metagga': 'MBJ', 

875 }, 

876 # vdW-DFs 

877 'vdw-df': { 

878 'gga': 'RE', 

879 'luse_vdw': True, 

880 'aggac': 0. 

881 }, 

882 'vdw-df-cx': { 

883 'gga': 'CX', 

884 'luse_vdw': True, 

885 'aggac': 0. 

886 }, 

887 'vdw-df-cx0p': { 

888 'gga': 'CX', 

889 'luse_vdw': True, 

890 'aggac': 0., 

891 'lhfcalc': True, 

892 'aexx': 0.2, 

893 'aggax': 0.8 

894 }, 

895 'optpbe-vdw': { 

896 'gga': 'OR', 

897 'luse_vdw': True, 

898 'aggac': 0.0 

899 }, 

900 'optb88-vdw': { 

901 'gga': 'BO', 

902 'luse_vdw': True, 

903 'aggac': 0.0, 

904 'param1': 1.1 / 6.0, 

905 'param2': 0.22 

906 }, 

907 'optb86b-vdw': { 

908 'gga': 'MK', 

909 'luse_vdw': True, 

910 'aggac': 0.0, 

911 'param1': 0.1234, 

912 'param2': 1.0 

913 }, 

914 'vdw-df2': { 

915 'gga': 'ML', 

916 'luse_vdw': True, 

917 'aggac': 0.0, 

918 'zab_vdw': -1.8867 

919 }, 

920 'rev-vdw-df2': { 

921 'gga': 'MK', 

922 'luse_vdw': True, 

923 'param1': 0.1234, 

924 'param2': 0.711357, 

925 'zab_vdw': -1.8867, 

926 'aggac': 0.0 

927 }, 

928 'beef-vdw': { 

929 'gga': 'BF', 

930 'luse_vdw': True, 

931 'zab_vdw': -1.8867 

932 }, 

933 # Hartree-Fock and hybrids 

934 'hf': { 

935 'lhfcalc': True, 

936 'aexx': 1.0, 

937 'aldac': 0.0, 

938 'aggac': 0.0 

939 }, 

940 'b3lyp': { 

941 'gga': 'B3', 

942 'lhfcalc': True, 

943 'aexx': 0.2, 

944 'aggax': 0.72, 

945 'aggac': 0.81, 

946 'aldac': 0.19 

947 }, 

948 'pbe0': { 

949 'gga': 'PE', 

950 'lhfcalc': True 

951 }, 

952 'hse03': { 

953 'gga': 'PE', 

954 'lhfcalc': True, 

955 'hfscreen': 0.3 

956 }, 

957 'hse06': { 

958 'gga': 'PE', 

959 'lhfcalc': True, 

960 'hfscreen': 0.2 

961 }, 

962 'hsesol': { 

963 'gga': 'PS', 

964 'lhfcalc': True, 

965 'hfscreen': 0.2 

966 }, 

967 # MN-VFM functionals 

968 'sogga': { 

969 'gga': 'SA' 

970 }, 

971 'sogga11': { 

972 'gga': 'S1' 

973 }, 

974 'sogga11-x': { 

975 'gga': 'SX', 

976 'lhfcalc': True, 

977 'aexx': 0.401 

978 }, 

979 'n12': { 

980 'gga': 'N2' 

981 }, 

982 'n12-sx': { 

983 'gga': 'NX', 

984 'lhfcalc': True, 

985 'lhfscreen': 0.2 

986 }, 

987 'mn12l': { 

988 'metagga': 'MN12L' 

989 }, 

990 'gam': { 

991 'gga': 'GA' 

992 }, 

993 'mn15l': { 

994 'metagga': 'MN15L' 

995 }, 

996 'hle17': { 

997 'metagga': 'HLE17' 

998 }, 

999 'revm06l': { 

1000 'metagga': 'revM06L' 

1001 }, 

1002 'm06sx': { 

1003 'metagga': 'M06SX', 

1004 'lhfcalc': True, 

1005 'hfscreen': 0.189, 

1006 'aexx': 0.335 

1007 } 

1008 } 

1009 

1010 # environment variable for PP paths 

1011 VASP_PP_PATH = 'VASP_PP_PATH' 

1012 

1013 def __init__(self, restart=None): 

1014 self.float_params = {} 

1015 self.exp_params = {} 

1016 self.string_params = {} 

1017 self.int_params = {} 

1018 self.bool_params = {} 

1019 self.list_bool_params = {} 

1020 self.list_int_params = {} 

1021 self.list_float_params = {} 

1022 self.special_params = {} 

1023 self.dict_params = {} 

1024 for key in float_keys: 

1025 self.float_params[key] = None 

1026 for key in exp_keys: 

1027 self.exp_params[key] = None 

1028 for key in string_keys: 

1029 self.string_params[key] = None 

1030 for key in int_keys: 

1031 self.int_params[key] = None 

1032 for key in bool_keys: 

1033 self.bool_params[key] = None 

1034 for key in list_bool_keys: 

1035 self.list_bool_params[key] = None 

1036 for key in list_int_keys: 

1037 self.list_int_params[key] = None 

1038 for key in list_float_keys: 

1039 self.list_float_params[key] = None 

1040 for key in special_keys: 

1041 self.special_params[key] = None 

1042 for key in dict_keys: 

1043 self.dict_params[key] = None 

1044 

1045 # Initialize internal dictionary of input parameters which are 

1046 # not regular VASP keys 

1047 self.input_params = { 

1048 'xc': None, # Exchange-correlation recipe (e.g. 'B3LYP') 

1049 'pp': None, # Pseudopotential file (e.g. 'PW91') 

1050 'setups': None, # Special setups (e.g pv, sv, ...) 

1051 'txt': '-', # Where to send information 

1052 'kpts': (1, 1, 1), # k-points 

1053 # Option to use gamma-sampling instead of Monkhorst-Pack: 

1054 'gamma': False, 

1055 # number of points between points in band structures: 

1056 'kpts_nintersections': None, 

1057 # Option to write explicit k-points in units 

1058 # of reciprocal lattice vectors: 

1059 'reciprocal': False, 

1060 # Switch to disable writing constraints to POSCAR 

1061 'ignore_constraints': False, 

1062 # Net charge for the whole system; determines nelect if not 0 

1063 'charge': None, 

1064 # Deprecated older parameter which works just like "charge" but 

1065 # with the sign flipped 

1066 'net_charge': None, 

1067 # Custom key-value pairs, written to INCAR with *no* type checking 

1068 'custom': {}, 

1069 } 

1070 

1071 def set_xc_params(self, xc): 

1072 """Set parameters corresponding to XC functional""" 

1073 xc = xc.lower() 

1074 if xc is None: 

1075 pass 

1076 elif xc not in self.xc_defaults: 

1077 xc_allowed = ', '.join(self.xc_defaults.keys()) 

1078 raise ValueError('{0} is not supported for xc! Supported xc values' 

1079 'are: {1}'.format(xc, xc_allowed)) 

1080 else: 

1081 # XC defaults to PBE pseudopotentials 

1082 if 'pp' not in self.xc_defaults[xc]: 

1083 self.set(pp='PBE') 

1084 self.set(**self.xc_defaults[xc]) 

1085 

1086 def set(self, **kwargs): 

1087 

1088 if (('ldauu' in kwargs) and ('ldaul' in kwargs) and ('ldauj' in kwargs) 

1089 and ('ldau_luj' in kwargs)): 

1090 raise NotImplementedError( 

1091 'You can either specify ldaul, ldauu, and ldauj OR ' 

1092 'ldau_luj. ldau_luj is not a VASP keyword. It is a ' 

1093 'dictionary that specifies L, U and J for each ' 

1094 'chemical species in the atoms object. ' 

1095 'For example for a water molecule:' 

1096 '''ldau_luj={'H':{'L':2, 'U':4.0, 'J':0.9}, 

1097 'O':{'L':2, 'U':4.0, 'J':0.9}}''') 

1098 

1099 if 'xc' in kwargs: 

1100 self.set_xc_params(kwargs['xc']) 

1101 for key in kwargs: 

1102 if key in self.float_params: 

1103 self.float_params[key] = kwargs[key] 

1104 elif key in self.exp_params: 

1105 self.exp_params[key] = kwargs[key] 

1106 elif key in self.string_params: 

1107 self.string_params[key] = kwargs[key] 

1108 elif key in self.int_params: 

1109 self.int_params[key] = kwargs[key] 

1110 elif key in self.bool_params: 

1111 self.bool_params[key] = kwargs[key] 

1112 elif key in self.list_bool_params: 

1113 self.list_bool_params[key] = kwargs[key] 

1114 elif key in self.list_int_params: 

1115 self.list_int_params[key] = kwargs[key] 

1116 elif key in self.list_float_params: 

1117 self.list_float_params[key] = kwargs[key] 

1118 elif key in self.special_params: 

1119 self.special_params[key] = kwargs[key] 

1120 elif key in self.dict_params: 

1121 self.dict_params[key] = kwargs[key] 

1122 elif key in self.input_params: 

1123 self.input_params[key] = kwargs[key] 

1124 else: 

1125 raise TypeError('Parameter not defined: ' + key) 

1126 

1127 def check_xc(self): 

1128 """Make sure the calculator has functional & pseudopotentials set up 

1129 

1130 If no XC combination, GGA functional or POTCAR type is specified, 

1131 default to PW91. Otherwise, try to guess the desired pseudopotentials. 

1132 """ 

1133 

1134 p = self.input_params 

1135 

1136 # There is no way to correctly guess the desired 

1137 # set of pseudopotentials without 'pp' being set. 

1138 # Usually, 'pp' will be set by 'xc'. 

1139 if 'pp' not in p or p['pp'] is None: 

1140 if self.string_params['gga'] is None: 

1141 p.update({'pp': 'lda'}) 

1142 elif self.string_params['gga'] == '91': 

1143 p.update({'pp': 'pw91'}) 

1144 elif self.string_params['gga'] == 'PE': 

1145 p.update({'pp': 'pbe'}) 

1146 else: 

1147 raise NotImplementedError( 

1148 "Unable to guess the desired set of pseudopotential" 

1149 "(POTCAR) files. Please do one of the following: \n" 

1150 "1. Use the 'xc' parameter to define your XC functional." 

1151 "These 'recipes' determine the pseudopotential file as " 

1152 "well as setting the INCAR parameters.\n" 

1153 "2. Use the 'gga' settings None (default), 'PE' or '91'; " 

1154 "these correspond to LDA, PBE and PW91 respectively.\n" 

1155 "3. Set the POTCAR explicitly with the 'pp' flag. The " 

1156 "value should be the name of a folder on the VASP_PP_PATH" 

1157 ", and the aliases 'LDA', 'PBE' and 'PW91' are also" 

1158 "accepted.\n") 

1159 

1160 if (p['xc'] is not None and p['xc'].lower() == 'lda' 

1161 and p['pp'].lower() != 'lda'): 

1162 warnings.warn("XC is set to LDA, but PP is set to " 

1163 "{0}. \nThis calculation is using the {0} " 

1164 "POTCAR set. \n Please check that this is " 

1165 "really what you intended!" 

1166 "\n".format(p['pp'].upper())) 

1167 

1168 def _make_sort( 

1169 self, atoms: ase.Atoms, special_setups: Sequence[int] = () 

1170 ) -> Tuple[List[int], List[int]]: 

1171 symbols, _ = count_symbols(atoms, exclude=special_setups) 

1172 

1173 # Create sorting list 

1174 srt = [] # type: List[int] 

1175 srt.extend(special_setups) 

1176 

1177 for symbol in symbols: 

1178 for m, atom in enumerate(atoms): 

1179 if m in special_setups: 

1180 continue 

1181 if atom.symbol == symbol: 

1182 srt.append(m) 

1183 # Create the resorting list 

1184 resrt = list(range(len(srt))) 

1185 for n in range(len(resrt)): 

1186 resrt[srt[n]] = n 

1187 return srt, resrt 

1188 

1189 def _build_pp_list(self, 

1190 atoms, 

1191 setups=None, 

1192 special_setups: Sequence[int] = ()): 

1193 """Build the pseudopotential lists""" 

1194 

1195 p = self.input_params 

1196 

1197 if setups is None: 

1198 setups, special_setups = self._get_setups() 

1199 

1200 symbols, _ = count_symbols(atoms, exclude=special_setups) 

1201 

1202 # Potpaw folders may be identified by an alias or full name 

1203 for pp_alias, pp_folder in (('lda', 'potpaw'), ('pw91', 'potpaw_GGA'), 

1204 ('pbe', 'potpaw_PBE')): 

1205 if p['pp'].lower() == pp_alias: 

1206 break 

1207 else: 

1208 pp_folder = p['pp'] 

1209 

1210 if self.VASP_PP_PATH in os.environ: 

1211 pppaths = os.environ[self.VASP_PP_PATH].split(':') 

1212 else: 

1213 pppaths = [] 

1214 ppp_list = [] 

1215 # Setting the pseudopotentials, first special setups and 

1216 # then according to symbols 

1217 for m in special_setups: 

1218 if m in setups: 

1219 special_setup_index = m 

1220 elif str(m) in setups: 

1221 special_setup_index = str(m) # type: ignore 

1222 else: 

1223 raise Exception("Having trouble with special setup index {0}." 

1224 " Please use an int.".format(m)) 

1225 potcar = join(pp_folder, setups[special_setup_index], 'POTCAR') 

1226 for path in pppaths: 

1227 filename = join(path, potcar) 

1228 

1229 if isfile(filename) or islink(filename): 

1230 ppp_list.append(filename) 

1231 break 

1232 elif isfile(filename + '.Z') or islink(filename + '.Z'): 

1233 ppp_list.append(filename + '.Z') 

1234 break 

1235 else: 

1236 symbol = atoms.symbols[m] 

1237 msg = """Looking for {}. 

1238 No pseudopotential for symbol{} with setup {} """.format( 

1239 potcar, symbol, setups[special_setup_index]) 

1240 raise RuntimeError(msg) 

1241 

1242 for symbol in symbols: 

1243 try: 

1244 potcar = join(pp_folder, symbol + setups[symbol], 'POTCAR') 

1245 except (TypeError, KeyError): 

1246 potcar = join(pp_folder, symbol, 'POTCAR') 

1247 for path in pppaths: 

1248 filename = join(path, potcar) 

1249 

1250 if isfile(filename) or islink(filename): 

1251 ppp_list.append(filename) 

1252 break 

1253 elif isfile(filename + '.Z') or islink(filename + '.Z'): 

1254 ppp_list.append(filename + '.Z') 

1255 break 

1256 else: 

1257 msg = ("""Looking for PP for {} 

1258 The pseudopotentials are expected to be in: 

1259 LDA: $VASP_PP_PATH/potpaw/ 

1260 PBE: $VASP_PP_PATH/potpaw_PBE/ 

1261 PW91: $VASP_PP_PATH/potpaw_GGA/ 

1262 

1263 No pseudopotential for {}!""".format(potcar, symbol)) 

1264 raise RuntimeError(msg) 

1265 return ppp_list 

1266 

1267 def _get_setups(self): 

1268 p = self.input_params 

1269 

1270 special_setups = [] 

1271 

1272 # Default setup lists are available: 'minimal', 'recommended' and 'GW' 

1273 # These may be provided as a string e.g.:: 

1274 # 

1275 # calc = Vasp(setups='recommended') 

1276 # 

1277 # or in a dict with other specifications e.g.:: 

1278 # 

1279 # calc = Vasp(setups={'base': 'minimal', 'Ca': '_sv', 2: 'O_s'}) 

1280 # 

1281 # Where other keys are either atom identities or indices, and the 

1282 # corresponding values are suffixes or the full name of the setup 

1283 # folder, respectively. 

1284 

1285 # Avoid mutating the module dictionary, so we use a copy instead 

1286 # Note, it is a nested dict, so a regular copy is not enough 

1287 setups_defaults = get_default_setups() 

1288 

1289 # Default to minimal basis 

1290 if p['setups'] is None: 

1291 p['setups'] = {'base': 'minimal'} 

1292 

1293 # String shortcuts are initialised to dict form 

1294 elif isinstance(p['setups'], str): 

1295 if p['setups'].lower() in setups_defaults.keys(): 

1296 p['setups'] = {'base': p['setups']} 

1297 

1298 # Dict form is then queried to add defaults from setups.py. 

1299 if 'base' in p['setups']: 

1300 setups = setups_defaults[p['setups']['base'].lower()] 

1301 else: 

1302 setups = {} 

1303 

1304 # Override defaults with user-defined setups 

1305 if p['setups'] is not None: 

1306 setups.update(p['setups']) 

1307 

1308 for m in setups: 

1309 try: 

1310 special_setups.append(int(m)) 

1311 except ValueError: 

1312 pass 

1313 return setups, special_setups 

1314 

1315 def initialize(self, atoms): 

1316 """Initialize a VASP calculation 

1317 

1318 Constructs the POTCAR file (does not actually write it). 

1319 User should specify the PATH 

1320 to the pseudopotentials in VASP_PP_PATH environment variable 

1321 

1322 The pseudopotentials are expected to be in: 

1323 LDA: $VASP_PP_PATH/potpaw/ 

1324 PBE: $VASP_PP_PATH/potpaw_PBE/ 

1325 PW91: $VASP_PP_PATH/potpaw_GGA/ 

1326 

1327 if your pseudopotentials are somewhere else, or named 

1328 differently you may make symlinks at the paths above that 

1329 point to the right place. Alternatively, you may pass the full 

1330 name of a folder on the VASP_PP_PATH to the 'pp' parameter. 

1331 """ 

1332 

1333 self.check_xc() 

1334 self.atoms = atoms 

1335 self.all_symbols = atoms.get_chemical_symbols() 

1336 self.natoms = len(atoms) 

1337 

1338 self.spinpol = (atoms.get_initial_magnetic_moments().any() 

1339 or self.int_params['ispin'] == 2) 

1340 

1341 setups, special_setups = self._get_setups() 

1342 

1343 # Determine the number of atoms of each atomic species 

1344 # sorted after atomic species 

1345 symbols, symbolcount = count_symbols(atoms, exclude=special_setups) 

1346 self.sort, self.resort = self._make_sort(atoms, 

1347 special_setups=special_setups) 

1348 

1349 self.atoms_sorted = atoms[self.sort] 

1350 

1351 # Check if the necessary POTCAR files exists and 

1352 # create a list of their paths. 

1353 atomtypes = atoms.get_chemical_symbols() 

1354 self.symbol_count = [] 

1355 for m in special_setups: 

1356 self.symbol_count.append([atomtypes[m], 1]) 

1357 for m in symbols: 

1358 self.symbol_count.append([m, symbolcount[m]]) 

1359 

1360 # create pseudopotential list 

1361 self.ppp_list = self._build_pp_list(atoms, 

1362 setups=setups, 

1363 special_setups=special_setups) 

1364 

1365 self.converged = None 

1366 self.setups_changed = None 

1367 

1368 def default_nelect_from_ppp(self): 

1369 """ Get default number of electrons from ppp_list and symbol_count 

1370 

1371 "Default" here means that the resulting cell would be neutral. 

1372 """ 

1373 symbol_valences = [] 

1374 for filename in self.ppp_list: 

1375 with open_potcar(filename=filename) as ppp_file: 

1376 r = read_potcar_numbers_of_electrons(ppp_file) 

1377 symbol_valences.extend(r) 

1378 assert len(self.symbol_count) == len(symbol_valences) 

1379 default_nelect = 0 

1380 for ((symbol1, count), 

1381 (symbol2, valence)) in zip(self.symbol_count, symbol_valences): 

1382 assert symbol1 == symbol2 

1383 default_nelect += count * valence 

1384 return default_nelect 

1385 

1386 def write_input(self, atoms, directory='./'): 

1387 from ase.io.vasp import write_vasp 

1388 write_vasp(join(directory, 'POSCAR'), 

1389 self.atoms_sorted, 

1390 symbol_count=self.symbol_count, 

1391 ignore_constraints=self.input_params['ignore_constraints']) 

1392 self.write_incar(atoms, directory=directory) 

1393 self.write_potcar(directory=directory) 

1394 self.write_kpoints(atoms=atoms, directory=directory) 

1395 self.write_sort_file(directory=directory) 

1396 self.copy_vdw_kernel(directory=directory) 

1397 

1398 def copy_vdw_kernel(self, directory='./'): 

1399 """Method to copy the vdw_kernel.bindat file. 

1400 Set ASE_VASP_VDW environment variable to the vdw_kernel.bindat 

1401 folder location. Checks if LUSE_VDW is enabled, and if no location 

1402 for the vdW kernel is specified, a warning is issued.""" 

1403 

1404 vdw_env = 'ASE_VASP_VDW' 

1405 kernel = 'vdw_kernel.bindat' 

1406 dst = os.path.join(directory, kernel) 

1407 

1408 # No need to copy the file again 

1409 if isfile(dst): 

1410 return 

1411 

1412 if self.bool_params['luse_vdw']: 

1413 src = None 

1414 if vdw_env in os.environ: 

1415 src = os.path.join(os.environ[vdw_env], kernel) 

1416 

1417 if not src or not isfile(src): 

1418 warnings.warn( 

1419 ('vdW has been enabled, however no' 

1420 ' location for the {} file' 

1421 ' has been specified.' 

1422 ' Set {} environment variable to' 

1423 ' copy the vdW kernel.').format(kernel, vdw_env)) 

1424 else: 

1425 shutil.copyfile(src, dst) 

1426 

1427 def clean(self): 

1428 """Method which cleans up after a calculation. 

1429 

1430 The default files generated by Vasp will be deleted IF this 

1431 method is called. 

1432 

1433 """ 

1434 files = [ 

1435 'CHG', 'CHGCAR', 'POSCAR', 'INCAR', 'CONTCAR', 'DOSCAR', 

1436 'EIGENVAL', 'IBZKPT', 'KPOINTS', 'OSZICAR', 'OUTCAR', 'PCDAT', 

1437 'POTCAR', 'vasprun.xml', 'WAVECAR', 'XDATCAR', 'PROCAR', 

1438 'ase-sort.dat', 'LOCPOT', 'AECCAR0', 'AECCAR1', 'AECCAR2' 

1439 ] 

1440 for f in files: 

1441 try: 

1442 os.remove(f) 

1443 except OSError: 

1444 pass 

1445 

1446 def write_incar(self, atoms, directory='./', **kwargs): 

1447 """Writes the INCAR file.""" 

1448 p = self.input_params 

1449 # jrk 1/23/2015 I added this flag because this function has 

1450 # two places where magmoms get written. There is some 

1451 # complication when restarting that often leads to magmom 

1452 # getting written twice. this flag prevents that issue. 

1453 magmom_written = False 

1454 incar = open(join(directory, 'INCAR'), 'w') 

1455 incar.write('INCAR created by Atomic Simulation Environment\n') 

1456 for key, val in self.float_params.items(): 

1457 if key == 'nelect': 

1458 charge = p.get('charge') 

1459 # Handle deprecated net_charge parameter (remove at some point) 

1460 net_charge = p.get('net_charge') 

1461 if net_charge is not None: 

1462 warnings.warn( 

1463 '`net_charge`, which is given in units of ' 

1464 'the *negative* elementary charge (i.e., the opposite ' 

1465 'of what one normally calls charge) has been ' 

1466 'deprecated in favor of `charge`, which is given in ' 

1467 'units of the positive elementary charge as usual', 

1468 category=FutureWarning) 

1469 if charge is not None and charge != -net_charge: 

1470 raise ValueError( 

1471 "can't give both net_charge and charge") 

1472 charge = -net_charge 

1473 # We need to determine the nelect resulting from a given net 

1474 # charge in any case if it's != 0, but if nelect is 

1475 # additionally given explicitly, then we need to determine it 

1476 # even for net charge of 0 to check for conflicts 

1477 if charge is not None and (charge != 0 or val is not None): 

1478 default_nelect = self.default_nelect_from_ppp() 

1479 nelect_from_charge = default_nelect - charge 

1480 if val is not None and val != nelect_from_charge: 

1481 raise ValueError('incompatible input parameters: ' 

1482 'nelect=%s, but charge=%s ' 

1483 '(neutral nelect is %s)' % 

1484 (val, charge, default_nelect)) 

1485 val = nelect_from_charge 

1486 if val is not None: 

1487 incar.write(' %s = %5.6f\n' % (key.upper(), val)) 

1488 for key, val in self.exp_params.items(): 

1489 if val is not None: 

1490 incar.write(' %s = %5.2e\n' % (key.upper(), val)) 

1491 for key, val in self.string_params.items(): 

1492 if val is not None: 

1493 incar.write(' %s = %s\n' % (key.upper(), val)) 

1494 for key, val in self.int_params.items(): 

1495 if val is not None: 

1496 incar.write(' %s = %d\n' % (key.upper(), val)) 

1497 if key == 'ichain' and val > 0: 

1498 incar.write(' IBRION = 3\n POTIM = 0.0\n') 

1499 for key, val in self.int_params.items(): 

1500 if key == 'iopt' and val is None: 

1501 print('WARNING: optimization is ' 

1502 'set to LFBGS (IOPT = 1)') 

1503 incar.write(' IOPT = 1\n') 

1504 for key, val in self.exp_params.items(): 

1505 if key == 'ediffg' and val is None: 

1506 RuntimeError('Please set EDIFFG < 0') 

1507 

1508 for key, val in self.list_bool_params.items(): 

1509 if val is None: 

1510 pass 

1511 else: 

1512 incar.write(' %s = ' % key.upper()) 

1513 [incar.write('%s ' % _to_vasp_bool(x)) for x in val] 

1514 incar.write('\n') 

1515 

1516 for key, val in self.list_int_params.items(): 

1517 if val is None: 

1518 pass 

1519 elif key == 'ldaul' and (self.dict_params['ldau_luj'] is not None): 

1520 pass 

1521 else: 

1522 incar.write(' %s = ' % key.upper()) 

1523 [incar.write('%d ' % x) for x in val] 

1524 incar.write('\n') 

1525 

1526 for key, val in self.list_float_params.items(): 

1527 if val is None: 

1528 pass 

1529 elif ((key in ('ldauu', 'ldauj')) 

1530 and (self.dict_params['ldau_luj'] is not None)): 

1531 pass 

1532 elif key == 'magmom': 

1533 if not len(val) == len(atoms): 

1534 msg = ('Expected length of magmom tag to be' 

1535 ' {}, i.e. 1 value per atom, but got {}').format( 

1536 len(atoms), len(val)) 

1537 raise ValueError(msg) 

1538 

1539 # Check if user remembered to specify ispin 

1540 # note: we do not overwrite ispin if ispin=1 

1541 if not self.int_params['ispin']: 

1542 self.spinpol = True 

1543 incar.write(' ispin = 2\n'.upper()) 

1544 

1545 incar.write(' %s = ' % key.upper()) 

1546 magmom_written = True 

1547 # Work out compact a*x b*y notation and write in this form 

1548 # Assume 1 magmom per atom, ordered as our atoms object 

1549 val = np.array(val) 

1550 val = val[self.sort] # Order in VASP format 

1551 

1552 # Compactify the magmom list to symbol order 

1553 lst = [[1, val[0]]] 

1554 for n in range(1, len(val)): 

1555 if val[n] == val[n - 1]: 

1556 lst[-1][0] += 1 

1557 else: 

1558 lst.append([1, val[n]]) 

1559 incar.write(' '.join( 

1560 ['{:d}*{:.4f}'.format(mom[0], mom[1]) for mom in lst])) 

1561 incar.write('\n') 

1562 else: 

1563 incar.write(' %s = ' % key.upper()) 

1564 [incar.write('%.4f ' % x) for x in val] 

1565 incar.write('\n') 

1566 

1567 for key, val in self.bool_params.items(): 

1568 if val is not None: 

1569 incar.write(' %s = ' % key.upper()) 

1570 if val: 

1571 incar.write('.TRUE.\n') 

1572 else: 

1573 incar.write('.FALSE.\n') 

1574 for key, val in self.special_params.items(): 

1575 if val is not None: 

1576 incar.write(' %s = ' % key.upper()) 

1577 if key == 'lreal': 

1578 if isinstance(val, str): 

1579 incar.write(val + '\n') 

1580 elif isinstance(val, bool): 

1581 if val: 

1582 incar.write('.TRUE.\n') 

1583 else: 

1584 incar.write('.FALSE.\n') 

1585 for key, val in self.dict_params.items(): 

1586 if val is not None: 

1587 if key == 'ldau_luj': 

1588 # User didn't turn on LDAU tag. 

1589 # Only turn on if ldau is unspecified 

1590 if self.bool_params['ldau'] is None: 

1591 self.bool_params['ldau'] = True 

1592 # At this point we have already parsed our bool params 

1593 incar.write(' LDAU = .TRUE.\n') 

1594 llist = ulist = jlist = '' 

1595 for symbol in self.symbol_count: 

1596 # default: No +U 

1597 luj = val.get(symbol[0], {'L': -1, 'U': 0.0, 'J': 0.0}) 

1598 llist += ' %i' % luj['L'] 

1599 ulist += ' %.3f' % luj['U'] 

1600 jlist += ' %.3f' % luj['J'] 

1601 incar.write(' LDAUL =%s\n' % llist) 

1602 incar.write(' LDAUU =%s\n' % ulist) 

1603 incar.write(' LDAUJ =%s\n' % jlist) 

1604 

1605 if (self.spinpol and not magmom_written 

1606 # We don't want to write magmoms if they are all 0. 

1607 # but we could still be doing a spinpol calculation 

1608 and atoms.get_initial_magnetic_moments().any()): 

1609 if not self.int_params['ispin']: 

1610 incar.write(' ispin = 2\n'.upper()) 

1611 # Write out initial magnetic moments 

1612 magmom = atoms.get_initial_magnetic_moments()[self.sort] 

1613 # unpack magmom array if three components specified 

1614 if magmom.ndim > 1: 

1615 magmom = [item for sublist in magmom for item in sublist] 

1616 list = [[1, magmom[0]]] 

1617 for n in range(1, len(magmom)): 

1618 if magmom[n] == magmom[n - 1]: 

1619 list[-1][0] += 1 

1620 else: 

1621 list.append([1, magmom[n]]) 

1622 incar.write(' magmom = '.upper()) 

1623 [incar.write('%i*%.4f ' % (mom[0], mom[1])) for mom in list] 

1624 incar.write('\n') 

1625 

1626 # Custom key-value pairs, which receive no formatting 

1627 # Use the comment "# <Custom ASE key>" to denote such 

1628 # a custom key-value pair, as we cannot otherwise 

1629 # reliably and easily identify such non-standard entries 

1630 custom_kv_pairs = p.get('custom') 

1631 for key, value in custom_kv_pairs.items(): 

1632 incar.write(' {} = {} # <Custom ASE key>\n'.format( 

1633 key.upper(), value)) 

1634 incar.close() 

1635 

1636 def write_kpoints(self, atoms=None, directory='./', **kwargs): 

1637 """Writes the KPOINTS file.""" 

1638 

1639 if atoms is None: 

1640 atoms = self.atoms 

1641 

1642 # Don't write anything if KSPACING is being used 

1643 if self.float_params['kspacing'] is not None: 

1644 if self.float_params['kspacing'] > 0: 

1645 return 

1646 else: 

1647 raise ValueError("KSPACING value {0} is not allowable. " 

1648 "Please use None or a positive number." 

1649 "".format(self.float_params['kspacing'])) 

1650 

1651 kpointstring = format_kpoints( 

1652 kpts=self.input_params['kpts'], 

1653 atoms=atoms, 

1654 reciprocal=self.input_params['reciprocal'], 

1655 gamma=self.input_params['gamma']) 

1656 with open(join(directory, 'KPOINTS'), 'w') as kpoints: 

1657 kpoints.write(kpointstring) 

1658 

1659 def write_potcar(self, suffix="", directory='./'): 

1660 """Writes the POTCAR file.""" 

1661 

1662 with open(join(directory, 'POTCAR' + suffix), 'w') as potfile: 

1663 for filename in self.ppp_list: 

1664 with open_potcar(filename=filename) as ppp_file: 

1665 for line in ppp_file: 

1666 potfile.write(line) 

1667 

1668 def write_sort_file(self, directory='./'): 

1669 """Writes a sortings file. 

1670 

1671 This file contains information about how the atoms are sorted in 

1672 the first column and how they should be resorted in the second 

1673 column. It is used for restart purposes to get sorting right 

1674 when reading in an old calculation to ASE.""" 

1675 

1676 with open(join(directory, 'ase-sort.dat'), 'w') as fd: 

1677 for n in range(len(self.sort)): 

1678 fd.write('%5i %5i \n' % (self.sort[n], self.resort[n])) 

1679 

1680 # The below functions are used to restart a calculation 

1681 

1682 def read_incar(self, filename): 

1683 """Method that imports settings from INCAR file. 

1684 

1685 Typically named INCAR.""" 

1686 

1687 self.spinpol = False 

1688 with open(filename, 'r') as fd: 

1689 lines = fd.readlines() 

1690 

1691 for line in lines: 

1692 try: 

1693 # Make multiplication, comments, and parameters easier to spot 

1694 line = line.replace("*", " * ") 

1695 line = line.replace("=", " = ") 

1696 line = line.replace("#", "# ") 

1697 data = line.split() 

1698 # Skip empty and commented lines. 

1699 if len(data) == 0: 

1700 continue 

1701 elif data[0][0] in ['#', '!']: 

1702 continue 

1703 key = data[0].lower() 

1704 if '<Custom ASE key>' in line: 

1705 # This key was added with custom key-value pair formatting. 

1706 # Unconditionally add it, no type checking 

1707 # Get value between "=" and the comment, e.g. 

1708 # key = 1 2 3 # <Custom ASE key> 

1709 # value should be '1 2 3' 

1710 

1711 # Split at first occurence of "=" 

1712 value = line.split('=', 1)[1] 

1713 # First "#" denotes beginning of comment 

1714 # Add everything before comment as a string to custom dict 

1715 value = value.split('#', 1)[0].strip() 

1716 self.input_params['custom'][key] = value 

1717 elif key in float_keys: 

1718 self.float_params[key] = float(data[2]) 

1719 elif key in exp_keys: 

1720 self.exp_params[key] = float(data[2]) 

1721 elif key in string_keys: 

1722 self.string_params[key] = str(data[2]) 

1723 elif key in int_keys: 

1724 if key == 'ispin': 

1725 # JRK added. not sure why we would want to leave ispin 

1726 # out 

1727 self.int_params[key] = int(data[2]) 

1728 if int(data[2]) == 2: 

1729 self.spinpol = True 

1730 else: 

1731 self.int_params[key] = int(data[2]) 

1732 elif key in bool_keys: 

1733 if 'true' in data[2].lower(): 

1734 self.bool_params[key] = True 

1735 elif 'false' in data[2].lower(): 

1736 self.bool_params[key] = False 

1737 

1738 elif key in list_bool_keys: 

1739 self.list_bool_params[key] = [ 

1740 _from_vasp_bool(x) 

1741 for x in _args_without_comment(data[2:]) 

1742 ] 

1743 

1744 elif key in list_int_keys: 

1745 self.list_int_params[key] = [ 

1746 int(x) for x in _args_without_comment(data[2:]) 

1747 ] 

1748 

1749 elif key in list_float_keys: 

1750 if key == 'magmom': 

1751 lst = [] 

1752 i = 2 

1753 while i < len(data): 

1754 if data[i] in ["#", "!"]: 

1755 break 

1756 if data[i] == "*": 

1757 b = lst.pop() 

1758 i += 1 

1759 for j in range(int(b)): 

1760 lst.append(float(data[i])) 

1761 else: 

1762 lst.append(float(data[i])) 

1763 i += 1 

1764 self.list_float_params['magmom'] = lst 

1765 lst = np.array(lst) 

1766 if self.atoms is not None: 

1767 self.atoms.set_initial_magnetic_moments( 

1768 lst[self.resort]) 

1769 else: 

1770 data = _args_without_comment(data) 

1771 self.list_float_params[key] = [ 

1772 float(x) for x in data[2:] 

1773 ] 

1774 # elif key in list_keys: 

1775 # list = [] 

1776 # if key in ('dipol', 'eint', 'ferwe', 'ferdo', 

1777 # 'ropt', 'rwigs', 

1778 # 'ldauu', 'ldaul', 'ldauj', 'langevin_gamma'): 

1779 # for a in data[2:]: 

1780 # if a in ["!", "#"]: 

1781 # break 

1782 # list.append(float(a)) 

1783 # elif key in ('iband', 'kpuse', 'random_seed'): 

1784 # for a in data[2:]: 

1785 # if a in ["!", "#"]: 

1786 # break 

1787 # list.append(int(a)) 

1788 # self.list_params[key] = list 

1789 # if key == 'magmom': 

1790 # list = [] 

1791 # i = 2 

1792 # while i < len(data): 

1793 # if data[i] in ["#", "!"]: 

1794 # break 

1795 # if data[i] == "*": 

1796 # b = list.pop() 

1797 # i += 1 

1798 # for j in range(int(b)): 

1799 # list.append(float(data[i])) 

1800 # else: 

1801 # list.append(float(data[i])) 

1802 # i += 1 

1803 # self.list_params['magmom'] = list 

1804 # list = np.array(list) 

1805 # if self.atoms is not None: 

1806 # self.atoms.set_initial_magnetic_moments( 

1807 # list[self.resort]) 

1808 elif key in special_keys: 

1809 if key == 'lreal': 

1810 if 'true' in data[2].lower(): 

1811 self.special_params[key] = True 

1812 elif 'false' in data[2].lower(): 

1813 self.special_params[key] = False 

1814 else: 

1815 self.special_params[key] = data[2] 

1816 except KeyError: 

1817 raise IOError('Keyword "%s" in INCAR is' 

1818 'not known by calculator.' % key) 

1819 except IndexError: 

1820 raise IOError('Value missing for keyword "%s".' % key) 

1821 

1822 def read_kpoints(self, filename): 

1823 """Read kpoints file, typically named KPOINTS.""" 

1824 # If we used VASP builtin kspacing, 

1825 if self.float_params['kspacing'] is not None: 

1826 # Don't update kpts array 

1827 return 

1828 

1829 with open(filename, 'r') as fd: 

1830 lines = fd.readlines() 

1831 

1832 ktype = lines[2].split()[0].lower()[0] 

1833 if ktype in ['g', 'm', 'a']: 

1834 if ktype == 'g': 

1835 self.set(gamma=True) 

1836 kpts = np.array([int(lines[3].split()[i]) for i in range(3)]) 

1837 elif ktype == 'a': 

1838 kpts = np.array([int(lines[3].split()[i]) for i in range(1)]) 

1839 elif ktype == 'm': 

1840 kpts = np.array([int(lines[3].split()[i]) for i in range(3)]) 

1841 else: 

1842 if ktype in ['c', 'k']: 

1843 self.set(reciprocal=False) 

1844 else: 

1845 self.set(reciprocal=True) 

1846 kpts = np.array( 

1847 [list(map(float, line.split())) for line in lines[3:]]) 

1848 self.set(kpts=kpts) 

1849 

1850 def read_potcar(self, filename): 

1851 """ Read the pseudopotential XC functional from POTCAR file. 

1852 """ 

1853 

1854 # Search for key 'LEXCH' in POTCAR 

1855 xc_flag = None 

1856 with open(filename, 'r') as fd: 

1857 for line in fd: 

1858 key = line.split()[0].upper() 

1859 if key == 'LEXCH': 

1860 xc_flag = line.split()[-1].upper() 

1861 break 

1862 

1863 if xc_flag is None: 

1864 raise ValueError('LEXCH flag not found in POTCAR file.') 

1865 

1866 # Values of parameter LEXCH and corresponding XC-functional 

1867 xc_dict = {'PE': 'PBE', '91': 'PW91', 'CA': 'LDA'} 

1868 

1869 if xc_flag not in xc_dict.keys(): 

1870 raise ValueError('Unknown xc-functional flag found in POTCAR,' 

1871 ' LEXCH=%s' % xc_flag) 

1872 

1873 self.input_params['pp'] = xc_dict[xc_flag] 

1874 

1875 def todict(self): 

1876 """Returns a dictionary of all parameters 

1877 that can be used to construct a new calculator object""" 

1878 dict_list = [ 

1879 'float_params', 'exp_params', 'string_params', 'int_params', 

1880 'bool_params', 'list_bool_params', 'list_int_params', 

1881 'list_float_params', 'special_params', 'dict_params', 

1882 'input_params' 

1883 ] 

1884 dct = {} 

1885 for item in dict_list: 

1886 dct.update(getattr(self, item)) 

1887 dct = {key: value for key, value in dct.items() if value is not None} 

1888 return dct 

1889 

1890 

1891def _args_without_comment(data, marks=['!', '#']): 

1892 """Check split arguments list for a comment, return data up to marker 

1893 

1894 INCAR reader splits list arguments on spaces and leaves comment markers as 

1895 individual items. This function returns only the data portion of the list. 

1896 

1897 """ 

1898 comment_locs = [data.index(mark) for mark in marks if mark in data] 

1899 if comment_locs == []: 

1900 return data 

1901 else: 

1902 return data[:min(comment_locs)] 

1903 

1904 

1905def _from_vasp_bool(x): 

1906 """Cast vasp boolean to Python bool 

1907 

1908 VASP files sometimes use T or F as shorthand for the preferred Boolean 

1909 notation .TRUE. or .FALSE. As capitalisation is pretty inconsistent in 

1910 practice, we allow all cases to be cast to a Python bool. 

1911 

1912 """ 

1913 assert isinstance(x, str) 

1914 if x.lower() == '.true.' or x.lower() == 't': 

1915 return True 

1916 elif x.lower() == '.false.' or x.lower() == 'f': 

1917 return False 

1918 else: 

1919 raise ValueError('Value "%s" not recognized as bool' % x) 

1920 

1921 

1922def _to_vasp_bool(x): 

1923 """Convert Python boolean to string for VASP input 

1924 

1925 In case the value was modified to a string already, appropriate strings 

1926 will also be accepted and cast to a standard .TRUE. / .FALSE. format. 

1927 

1928 """ 

1929 if isinstance(x, str): 

1930 if x.lower() in ('.true.', 't'): 

1931 x = True 

1932 elif x.lower() in ('.false.', 'f'): 

1933 x = False 

1934 else: 

1935 raise ValueError('"%s" not recognised as VASP Boolean') 

1936 assert isinstance(x, bool) 

1937 if x: 

1938 return '.TRUE.' 

1939 else: 

1940 return '.FALSE.' 

1941 

1942 

1943def open_potcar(filename): 

1944 """ Open POTCAR file with transparent decompression if it's an archive (.Z) 

1945 """ 

1946 import gzip 

1947 if filename.endswith('R'): 

1948 return open(filename, 'r') 

1949 elif filename.endswith('.Z'): 

1950 return gzip.open(filename) 

1951 else: 

1952 raise ValueError('Invalid POTCAR filename: "%s"' % filename) 

1953 

1954 

1955def read_potcar_numbers_of_electrons(file_obj): 

1956 """ Read list of tuples (atomic symbol, number of valence electrons) 

1957 for each atomtype from a POTCAR file.""" 

1958 nelect = [] 

1959 lines = file_obj.readlines() 

1960 for n, line in enumerate(lines): 

1961 if 'TITEL' in line: 

1962 symbol = line.split('=')[1].split()[1].split('_')[0].strip() 

1963 valence = float( 

1964 lines[n + 4].split(';')[1].split('=')[1].split()[0].strip()) 

1965 nelect.append((symbol, valence)) 

1966 return nelect 

1967 

1968 

1969def count_symbols(atoms, exclude=()): 

1970 """Count symbols in atoms object, excluding a set of indices 

1971 

1972 Parameters: 

1973 atoms: Atoms object to be grouped 

1974 exclude: List of indices to be excluded from the counting 

1975 

1976 Returns: 

1977 Tuple of (symbols, symbolcount) 

1978 symbols: The unique symbols in the included list 

1979 symbolscount: Count of symbols in the included list 

1980 

1981 Example: 

1982 

1983 >>> from ase.build import bulk 

1984 >>> atoms = bulk('NaCl', crystalstructure='rocksalt', a=4.1, cubic=True) 

1985 >>> count_symbols(atoms) 

1986 (['Na', 'Cl'], {'Na': 4, 'Cl': 4}) 

1987 >>> count_symbols(atoms, exclude=(1, 2, 3)) 

1988 (['Na', 'Cl'], {'Na': 3, 'Cl': 2}) 

1989 """ 

1990 symbols = [] 

1991 symbolcount = {} 

1992 for m, symbol in enumerate(atoms.symbols): 

1993 if m in exclude: 

1994 continue 

1995 if symbol not in symbols: 

1996 symbols.append(symbol) 

1997 symbolcount[symbol] = 1 

1998 else: 

1999 symbolcount[symbol] += 1 

2000 return symbols, symbolcount