1# coding: utf-8 2"""Factory functions for Abinit input files """ 3import numpy as np 4import pymatgen.io.abinit.abiobjects as aobj 5import abipy.abio.input_tags as atags 6 7from enum import Enum 8from collections import namedtuple 9from monty.collections import AttrDict 10from monty.string import is_string 11from monty.json import jsanitize, MontyDecoder, MSONable 12from pymatgen.util.serialization import pmg_serialize 13from abipy.core.structure import Structure 14from abipy.abio.inputs import AbinitInput, MultiDataset 15 16 17__all__ = [ 18 "gs_input", 19 "ebands_input", 20 "phonons_from_gsinput", 21 "g0w0_with_ppmodel_inputs", 22 "g0w0_convergence_inputs", 23 "bse_with_mdf_inputs", 24 "ion_ioncell_relax_input", 25 "ion_ioncell_relax_and_ebands_input", 26 "scf_phonons_inputs", 27 "piezo_elastic_inputs_from_gsinput", 28 "scf_piezo_elastic_inputs", 29 "scf_for_phonons", 30 "dte_from_gsinput", 31 "dfpt_from_gsinput", 32 "minimal_scf_input" 33] 34 35 36# Name of the (default) tolerance used by the runlevels. 37_runl2tolname = { 38 "scf": 'tolvrs', 39 "nscf": 'tolwfr', 40 "dfpt": 'toldfe', # ? 41 "screening": 'toldfe', # dummy 42 "sigma": 'toldfe', # dummy 43 "bse": 'toldfe', # ? 44 "relax": 'tolrff', 45} 46 47# Tolerances for the different levels of accuracy. 48T = namedtuple('Tolerance', "low normal high") 49_tolerances = { 50 "toldfe": T(1.e-7, 1.e-8, 1.e-9), 51 "tolvrs": T(1.e-7, 1.e-8, 1.e-9), 52 "tolwfr": T(1.e-15, 1.e-17, 1.e-19), 53 "tolrff": T(0.04, 0.02, 0.01)} 54del T 55 56 57# Default values used if user does not specify them 58_DEFAULTS = dict( 59 kppa=1000, 60) 61 62 63class ShiftMode(Enum): 64 """ 65 Class defining the mode to be used for the shifts. 66 G: Gamma centered 67 M: Monkhorst-Pack ((0.5, 0.5, 0.5)) 68 S: Symmetric. Respects the chksymbreak with multiple shifts 69 O: OneSymmetric. Respects the chksymbreak with a single shift (as in 'S' if a single shift is given, gamma 70 centered otherwise. 71 """ 72 GammaCentered = 'G' 73 MonkhorstPack = 'M' 74 Symmetric = 'S' 75 OneSymmetric = 'O' 76 77 @classmethod 78 def from_object(cls, obj): 79 """ 80 Returns an instance of ShiftMode based on the type of object passed. Converts strings to ShiftMode depending 81 on the iniital letter of the string. G for GammaCenterd, M for MonkhorstPack, S for Symmetric, O for OneSymmetric. 82 Case insensitive. 83 """ 84 if isinstance(obj, cls): 85 return obj 86 elif is_string(obj): 87 return cls(obj[0].upper()) 88 else: 89 raise TypeError('The object provided is not handled: type %s' % type(obj)) 90 91 92def _stopping_criterion(runlevel, accuracy): 93 """Return the stopping criterion for this runlevel with the given accuracy.""" 94 tolname = _runl2tolname[runlevel] 95 return {tolname: getattr(_tolerances[tolname], accuracy)} 96 97 98def _find_ecut_pawecutdg(ecut, pawecutdg, pseudos, accuracy): 99 """Return a |AttrDict| with the value of ``ecut`` and ``pawecutdg``.""" 100 # Get ecut and pawecutdg from the pseudo hints. 101 if ecut is None or (pawecutdg is None and any(p.ispaw for p in pseudos)): 102 has_hints = all(p.has_hints for p in pseudos) 103 104 if ecut is None: 105 if has_hints: 106 ecut = max(p.hint_for_accuracy(accuracy).ecut for p in pseudos) 107 else: 108 raise AbinitInput.Error("ecut is None but pseudos do not provide hints for ecut") 109 110 if pawecutdg is None and any(p.ispaw for p in pseudos): 111 if has_hints: 112 pawecutdg = max(p.hint_for_accuracy(accuracy).pawecutdg for p in pseudos) 113 else: 114 raise RuntimeError("pawecutdg is None but pseudos do not provide hints") 115 116 return AttrDict(ecut=ecut, pawecutdg=pawecutdg) 117 118 119def _find_scf_nband(structure, pseudos, electrons, spinat=None): 120 """Find the value of ``nband``.""" 121 if electrons.nband is not None: return electrons.nband 122 123 nsppol, smearing = electrons.nsppol, electrons.smearing 124 125 # Number of valence electrons including possible extra charge 126 nval = structure.num_valence_electrons(pseudos) 127 nval -= electrons.charge 128 129 # First guess (semiconductors) 130 nband = nval // 2 131 132 # TODO: Find better algorithm 133 # If nband is too small we may kill the job, increase nband and restart 134 # but this change could cause problems in the other steps of the calculation 135 # if the change is not propagated e.g. phonons in metals. 136 if smearing: 137 # metallic occupation 138 nband = max(np.ceil(nband * 1.2), nband + 10) 139 else: 140 nband = max(np.ceil(nband * 1.1), nband + 4) 141 142 # Increase number of bands based on the starting magnetization 143 if nsppol == 2 and spinat is not None: 144 nband += np.ceil(max(np.sum(spinat, axis=0)) / 2.) 145 146 # Force even nband (easier to divide among procs, mandatory if nspinor == 2) 147 nband += nband % 2 148 return int(nband) 149 150 151def _get_shifts(shift_mode, structure): 152 """ 153 Gives the shifts based on the selected shift mode and on the symmetry of the structure. 154 G: Gamma centered 155 M: Monkhorst-Pack ((0.5, 0.5, 0.5)) 156 S: Symmetric. Respects the chksymbreak with multiple shifts 157 O: OneSymmetric. Respects the chksymbreak with a single shift (as in 'S' if a single shift is given, gamma 158 centered otherwise. 159 160 Note: for some cases (e.g. body centered tetragonal), both the Symmetric and OneSymmetric may fail to satisfy the 161 ``chksymbreak`` condition (Abinit input variable). 162 """ 163 if shift_mode == ShiftMode.GammaCentered: 164 return ((0, 0, 0)) 165 elif shift_mode == ShiftMode.MonkhorstPack: 166 return ((0.5, 0.5, 0.5)) 167 elif shift_mode == ShiftMode.Symmetric: 168 structure = Structure.from_sites(structure) 169 return structure.calc_shiftk() 170 elif shift_mode == ShiftMode.OneSymmetric: 171 structure = Structure.from_sites(structure) 172 shifts = structure.calc_shiftk() 173 if len(shifts) == 1: 174 return shifts 175 else: 176 return ((0, 0, 0)) 177 else: 178 raise ValueError("invalid shift_mode: `%s`" % str(shift_mode)) 179 180 181def gs_input(structure, pseudos, 182 kppa=None, ecut=None, pawecutdg=None, scf_nband=None, accuracy="normal", spin_mode="polarized", 183 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None): 184 """ 185 Returns a |AbinitInput| for ground-state calculation. 186 187 Args: 188 structure: |Structure| object. 189 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 190 kppa: Defines the sampling used for the SCF run. Defaults to 1000 if not given. 191 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 192 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized from the pseudos 193 according to accuracy) 194 scf_nband: Number of bands for SCF run. If scf_nband is None, nband is automatically initialized 195 from the list of pseudos, the structure and the smearing option. 196 accuracy: Accuracy of the calculation. 197 spin_mode: Spin polarization. 198 smearing: Smearing technique. 199 charge: Electronic charge added to the unit cell. 200 scf_algorithm: Algorithm used for solving of the SCF cycle. 201 """ 202 multi = ebands_input(structure, pseudos, 203 kppa=kppa, ndivsm=0, 204 ecut=ecut, pawecutdg=pawecutdg, scf_nband=scf_nband, accuracy=accuracy, spin_mode=spin_mode, 205 smearing=smearing, charge=charge, scf_algorithm=scf_algorithm) 206 207 return multi[0] 208 209 210def ebands_input(structure, pseudos, 211 kppa=None, nscf_nband=None, ndivsm=15, 212 ecut=None, pawecutdg=None, scf_nband=None, accuracy="normal", spin_mode="polarized", 213 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None, dos_kppa=None): 214 """ 215 Returns a |MultiDataset| object for band structure calculations. 216 217 Args: 218 structure: |Structure| object. 219 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 220 kppa: Defines the sampling used for the SCF run. Defaults to 1000 if not given. 221 nscf_nband: Number of bands included in the NSCF run. Set to scf_nband + 10 if None. 222 ndivsm: Number of divisions used to sample the smallest segment of the k-path. 223 if 0, only the GS input is returned in multi[0]. 224 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 225 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized from the pseudos 226 according to accuracy) 227 scf_nband: Number of bands for SCF run. If scf_nband is None, nband is automatically initialized 228 from the list of pseudos, the structure and the smearing option. 229 accuracy: Accuracy of the calculation. 230 spin_mode: Spin polarization. 231 smearing: Smearing technique. 232 charge: Electronic charge added to the unit cell. 233 scf_algorithm: Algorithm used for solving of the SCF cycle. 234 dos_kppa: Scalar or List of integers with the number of k-points per atom 235 to be used for the computation of the DOS (None if DOS is not wanted). 236 """ 237 structure = Structure.as_structure(structure) 238 239 if dos_kppa is not None and not isinstance(dos_kppa, (list, tuple)): 240 dos_kppa = [dos_kppa] 241 242 multi = MultiDataset(structure, pseudos, ndtset=2 if dos_kppa is None else 2 + len(dos_kppa)) 243 244 # Set the cutoff energies. 245 multi.set_vars(_find_ecut_pawecutdg(ecut, pawecutdg, multi.pseudos, accuracy)) 246 247 # SCF calculation. 248 kppa = _DEFAULTS.get("kppa") if kppa is None else kppa 249 scf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0) 250 scf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 251 charge=charge, nband=scf_nband, fband=None) 252 253 if spin_mode == "polarized": 254 multi[0].set_autospinat() 255 256 if scf_electrons.nband is None: 257 scf_electrons.nband = _find_scf_nband(structure, multi.pseudos, scf_electrons, multi[0].get('spinat', None)) 258 259 multi[0].set_vars(scf_ksampling.to_abivars()) 260 multi[0].set_vars(scf_electrons.to_abivars()) 261 multi[0].set_vars(_stopping_criterion("scf", accuracy)) 262 if ndivsm == 0: return multi 263 264 # Band structure calculation. 265 nscf_ksampling = aobj.KSampling.path_from_structure(ndivsm, structure) 266 nscf_nband = scf_electrons.nband + 10 if nscf_nband is None else nscf_nband 267 nscf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm={"iscf": -2}, 268 charge=charge, nband=nscf_nband, fband=None) 269 270 multi[1].set_vars(nscf_ksampling.to_abivars()) 271 multi[1].set_vars(nscf_electrons.to_abivars()) 272 multi[1].set_vars(_stopping_criterion("nscf", accuracy)) 273 274 # DOS calculation with different values of kppa. 275 if dos_kppa is not None: 276 for i, kppa in enumerate(dos_kppa): 277 dos_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0) 278 #dos_ksampling = aobj.KSampling.monkhorst(dos_ngkpt, shiftk=dos_shiftk, chksymbreak=0) 279 dos_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm={"iscf": -2}, 280 charge=charge, nband=nscf_nband) 281 dt = 2 + i 282 multi[dt].set_vars(dos_ksampling.to_abivars()) 283 multi[dt].set_vars(dos_electrons.to_abivars()) 284 multi[dt].set_vars(_stopping_criterion("nscf", accuracy)) 285 286 return multi 287 288 289def ion_ioncell_relax_input(structure, pseudos, 290 kppa=None, nband=None, 291 ecut=None, pawecutdg=None, accuracy="normal", spin_mode="polarized", 292 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None, shift_mode='Monkhorst-pack'): 293 """ 294 Returns a |MultiDataset| for a structural relaxation. The first dataset optmizes the 295 atomic positions at fixed unit cell. The second datasets optimizes both ions and unit cell parameters. 296 297 Args: 298 structure: |Structure| object. 299 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 300 kppa: Defines the sampling used for the Brillouin zone. 301 nband: Number of bands included in the SCF run. 302 accuracy: Accuracy of the calculation. 303 spin_mode: Spin polarization. 304 smearing: Smearing technique. 305 charge: Electronic charge added to the unit cell. 306 scf_algorithm: Algorithm used for the solution of the SCF cycle. 307 """ 308 structure = Structure.as_structure(structure) 309 multi = MultiDataset(structure, pseudos, ndtset=2) 310 311 # Set the cutoff energies. 312 multi.set_vars(_find_ecut_pawecutdg(ecut, pawecutdg, multi.pseudos, accuracy)) 313 314 kppa = _DEFAULTS.get("kppa") if kppa is None else kppa 315 316 shift_mode = ShiftMode.from_object(shift_mode) 317 shifts = _get_shifts(shift_mode, structure) 318 ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=shifts) 319 electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 320 charge=charge, nband=nband, fband=None) 321 322 if spin_mode == "polarized": 323 spinat_dict = multi[0].set_autospinat() 324 multi[1].set_vars(spinat_dict) 325 326 if electrons.nband is None: 327 electrons.nband = _find_scf_nband(structure, multi.pseudos, electrons, multi[0].get('spinat', None)) 328 329 ion_relax = aobj.RelaxationMethod.atoms_only(atoms_constraints=None) 330 ioncell_relax = aobj.RelaxationMethod.atoms_and_cell(atoms_constraints=None) 331 332 multi.set_vars(electrons.to_abivars()) 333 multi.set_vars(ksampling.to_abivars()) 334 335 multi[0].set_vars(ion_relax.to_abivars()) 336 multi[0].set_vars(_stopping_criterion("relax", accuracy)) 337 338 multi[1].set_vars(ioncell_relax.to_abivars()) 339 multi[1].set_vars(_stopping_criterion("relax", accuracy)) 340 341 return multi 342 343 344def ion_ioncell_relax_and_ebands_input(structure, pseudos, 345 kppa=None, nband=None, 346 ecut=None, pawecutdg=None, accuracy="normal", spin_mode="polarized", 347 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None): 348 """ 349 Returns a |MultiDataset| for a structural relaxation followed by a band structure run. 350 The first dataset optimizes the atomic positions at fixed unit cell. 351 The second datasets optimizes both ions and unit cell parameters. 352 The other datasets perform a band structure calculation. 353 354 .. warning:: 355 356 Client code is responsible for propagating the relaxed structure obtained with the 357 second dataset to the inputs used for the band structure calculation. 358 359 Args: 360 structure: |Structure| object. 361 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 362 kppa: Defines the sampling used for the Brillouin zone. 363 nband: Number of bands included in the SCF run. 364 accuracy: Accuracy of the calculation. 365 spin_mode: Spin polarization. 366 smearing: Smearing technique. 367 charge: Electronic charge added to the unit cell. 368 scf_algorithm: Algorithm used for solving of the SCF cycle. 369 370 Returns: |MultiDataset| object 371 """ 372 structure = Structure.as_structure(structure) 373 374 relax_multi = ion_ioncell_relax_input(structure, pseudos, 375 kppa=kppa, nband=nband, 376 ecut=ecut, pawecutdg=pawecutdg, accuracy=accuracy, spin_mode=spin_mode, 377 smearing=smearing, charge=charge, scf_algorithm=scf_algorithm) 378 379 ebands_multi = ebands_input(structure, pseudos, 380 kppa=kppa, nscf_nband=None, ndivsm=15, 381 ecut=ecut, pawecutdg=pawecutdg, scf_nband=None, accuracy=accuracy, spin_mode=spin_mode, 382 smearing=smearing, charge=charge, scf_algorithm=scf_algorithm, dos_kppa=None) 383 384 return relax_multi + ebands_multi 385 386 387def g0w0_with_ppmodel_inputs(structure, pseudos, 388 kppa, nscf_nband, ecuteps, ecutsigx, 389 ecut=None, pawecutdg=None, shifts=(0.0, 0.0, 0.0), 390 accuracy="normal", spin_mode="polarized", smearing="fermi_dirac:0.1 eV", 391 ppmodel="godby", charge=0.0, scf_algorithm=None, inclvkb=2, scr_nband=None, 392 sigma_nband=None, gw_qprange=1): 393 """ 394 Returns a |MultiDataset| object that performs G0W0 calculations with the plasmon pole approximation. 395 396 Args: 397 structure: |Structure| object. 398 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 399 kppa: Defines the sampling used for the SCF run. 400 nscf_nband: Number of bands included in the NSCF run. 401 ecuteps: Cutoff energy [Ha] for the screening matrix. 402 ecutsigx: Cutoff energy [Ha] for the exchange part of the self-energy. 403 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 404 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized 405 from the pseudos according to accuracy) 406 shifts: Shifts for k-mesh. 407 accuracy: Accuracy of the calculation. 408 spin_mode: Spin polarization. 409 smearing: Smearing technique. 410 ppmodel: Plasmonpole technique. 411 charge: Electronic charge added to the unit cell. 412 scf_algorithm: Algorithm used for solving of the SCF cycle. 413 inclvkb: Treatment of the dipole matrix elements (see abinit variable). 414 scr_nband: Number of bands used to compute the screening (default is nscf_nband) 415 sigma_nband: Number of bands used to compute the self-energy (default is nscf_nband) 416 gw_qprange: Option for the automatic selection of k-points and bands for GW corrections. 417 See Abinit docs for more detail. The default value makes the code compute the 418 QP energies for all the point in the IBZ and one band above and one band below the Fermi level. 419 420 .. versionchanged: 0.3 421 422 The default value of ``shifts`` changed in v0.3 from (0.5, 0.5, 0.5) to (0.0, 0.0, 0.0). 423 """ 424 425 structure = Structure.as_structure(structure) 426 multi = MultiDataset(structure, pseudos, ndtset=4) 427 428 # Set the cutoff energies. 429 multi.set_vars(_find_ecut_pawecutdg(ecut, pawecutdg, multi.pseudos, accuracy)) 430 431 scf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=shifts) 432 scf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 433 charge=charge, nband=None, fband=None) 434 435 if scf_electrons.nband is None: 436 scf_electrons.nband = _find_scf_nband(structure, multi.pseudos, scf_electrons) 437 438 multi[0].set_vars(scf_ksampling.to_abivars()) 439 multi[0].set_vars(scf_electrons.to_abivars()) 440 multi[0].set_vars(_stopping_criterion("scf", accuracy)) 441 442 nscf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=shifts) 443 nscf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm={"iscf": -2}, 444 charge=charge, nband=nscf_nband, fband=None) 445 446 multi[1].set_vars(nscf_ksampling.to_abivars()) 447 multi[1].set_vars(nscf_electrons.to_abivars()) 448 multi[1].set_vars(_stopping_criterion("nscf", accuracy)) 449 # nbdbuf 450 451 # Screening. 452 if scr_nband is None: scr_nband = nscf_nband 453 screening = aobj.Screening(ecuteps, scr_nband, w_type="RPA", sc_mode="one_shot", 454 hilbert=None, ecutwfn=None, inclvkb=inclvkb) 455 456 multi[2].set_vars(nscf_ksampling.to_abivars()) 457 multi[2].set_vars(nscf_electrons.to_abivars()) 458 multi[2].set_vars(screening.to_abivars()) 459 multi[2].set_vars(_stopping_criterion("screening", accuracy)) # Dummy 460 461 # Sigma. 462 if sigma_nband is None: sigma_nband = nscf_nband 463 self_energy = aobj.SelfEnergy("gw", "one_shot", sigma_nband, ecutsigx, screening, 464 gw_qprange=gw_qprange, ppmodel=ppmodel) 465 466 multi[3].set_vars(nscf_ksampling.to_abivars()) 467 multi[3].set_vars(nscf_electrons.to_abivars()) 468 multi[3].set_vars(self_energy.to_abivars()) 469 multi[3].set_vars(_stopping_criterion("sigma", accuracy)) # Dummy 470 471 # TODO: Cannot use istwfk != 1. 472 multi.set_vars(istwfk="*1") 473 474 return multi 475 476 477def g0w0_convergence_inputs(structure, pseudos, kppa, nscf_nband, ecuteps, ecutsigx, scf_nband, ecut, 478 accuracy="normal", spin_mode="polarized", smearing="fermi_dirac:0.1 eV", 479 response_models=None, charge=0.0, scf_algorithm=None, inclvkb=2, 480 gw_qprange=1, gamma=True, nksmall=None, extra_abivars=None): 481 """ 482 Returns a |MultiDataset| object to generate a G0W0 work for the given the material. 483 See also :cite:`Setten2017`. 484 485 Args: 486 structure: |Structure| object 487 pseudos: List of |Pseudo| objects. 488 kppa: k points per reciprocal atom. 489 scf_nband: number of scf bands 490 ecut: ecut for all calcs that that are not ecut convergence cals at scf level 491 scf_ Defines the sampling used for the SCF run. 492 nscf_nband: a list of number of bands included in the screening and sigmaruns. 493 The NSCF run will be done with the maximum. 494 ecuteps: list of Cutoff energy [Ha] for the screening matrix. 495 ecutsigx: Cutoff energy [Ha] for the exchange part of the self-energy. 496 accuracy: Accuracy of the calculation. 497 spin_mode: Spin polarization. 498 smearing: Smearing technique. 499 charge: Electronic charge added to the unit cell. 500 scf_algorithm: Algorithm used for solving of the SCF cycle. 501 inclvkb: Treatment of the dipole matrix elements (see abinit variable). 502 response_models: List of response models 503 gw_qprange: selectpr for the qpoint mesh 504 gamma: is true a gamma centered mesh is enforced 505 nksmall: Kpoint division for additional band and dos calculations 506 extra_abivars: Dictionary with extra variables passed to ABINIT for all tasks. 507 508 extra abivars that are provided with _s appended will be take as a list of values to be tested a scf level 509 """ 510 if extra_abivars is None: 511 extra_abivars = {} 512 513 if response_models is None: 514 response_models = ["godby"] 515 516 scf_diffs = [] 517 518 keys = list(extra_abivars.keys()) 519 #for k in extra_abivars.keys(): 520 for k in keys: 521 if k[-2:] == '_s': 522 var = k[:len(k)-2] 523 values = extra_abivars.pop(k) 524 # to_add.update({k: values[-1]}) 525 for value in values: 526 diff_abivars = dict() 527 diff_abivars[var] = value 528 if pseudos.allpaw and var == 'ecut': 529 diff_abivars['pawecutdg'] = diff_abivars['ecut'] * 2 530 scf_diffs.append(diff_abivars) 531 532 extra_abivars_all = dict( 533 ecut=ecut, 534 paral_kgb=1, 535 istwfk="*1", 536 timopt=-1, 537 nbdbuf=8, 538 ) 539 540 extra_abivars_all.update(extra_abivars) 541 542 if pseudos.allpaw: 543 extra_abivars_all['pawecutdg'] = extra_abivars_all['ecut'] * 2 544 545 extra_abivars_gw = dict( 546 inclvkb=2, 547 symsigma=1, 548 gwpara=2, 549 gwmem='10', 550 prtsuscep=0 551 ) 552 553 # all these too many options are for development only the current idea for the final version is 554 #if gamma: 555 # scf_ksampling = aobj.KSampling.automatic_density(structure=structure, kppa=10000, chksymbreak=0, shifts=(0, 0, 0)) 556 # nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(2, 2, 2)) 557 # if kppa <= 13: 558 # nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(scf_kppa, scf_kppa, scf_kppa)) 559 # else: 560 # nscf_ksampling = aobj.KSampling.automatic_density(structure, scf_kppa, chksymbreak=0, shifts=(0, 0, 0)) 561 #else: 562 # scf_ksampling = aobj.KSampling.automatic_density(structure, scf_kppa, chksymbreak=0) 563 # nscf_ksampling = aobj.KSampling.automatic_density(structure, scf_kppa, chksymbreak=0) 564 565 if gamma: 566 if kppa == 1: 567 scf_ksampling = aobj.KSampling.gamma_centered(kpts=(1, 1, 1)) 568 nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(1, 1, 1)) 569 elif kppa == 2: 570 scf_ksampling = aobj.KSampling.gamma_centered(kpts=(2, 2, 2)) 571 nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(2, 2, 2)) 572 elif kppa < 0: 573 scf_ksampling = aobj.KSampling.gamma_centered(kpts=(-kppa, -kppa, -kppa)) 574 nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(2, 2, 2)) 575 elif kppa <= 13: 576 scf_ksampling = aobj.KSampling.gamma_centered(kpts=(kppa, kppa, kppa)) 577 nscf_ksampling = aobj.KSampling.gamma_centered(kpts=(kppa, kppa, kppa)) 578 else: 579 scf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=(0, 0, 0)) 580 nscf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=(0, 0, 0)) 581 else: 582 # this is the original behaviour before the development of the gwwrapper 583 scf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0) 584 nscf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0) 585 586 scf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 587 charge=charge, nband=scf_nband, fband=None) 588 nscf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm={"iscf": -2}, 589 charge=charge, nband=max(nscf_nband), fband=None) 590 591 multi_scf = MultiDataset(structure, pseudos, ndtset=max(1, len(scf_diffs))) 592 593 multi_scf.set_vars(scf_ksampling.to_abivars()) 594 multi_scf.set_vars(scf_electrons.to_abivars()) 595 multi_scf.set_vars(extra_abivars_all) 596 multi_scf.set_vars(_stopping_criterion(runlevel="scf", accuracy=accuracy)) 597 multi_scf.set_vars(extra_abivars) 598 599 for variables, abinput in zip(scf_diffs, multi_scf): 600 abinput.set_vars(variables) 601 602 scf_inputs = multi_scf.split_datasets() 603 604 # create nscf inputs 605 ndtset = 3 if nksmall is not None else 1 606 nscf_multi = MultiDataset(structure=structure, pseudos=pseudos, ndtset=ndtset) 607 608 nscf_multi.set_vars(nscf_electrons.to_abivars()) 609 nscf_multi.set_vars(extra_abivars_all) 610 nscf_multi.set_vars(_stopping_criterion(runlevel="nscf", accuracy=accuracy)) 611 612 nscf_multi[-1].set_vars(nscf_ksampling.to_abivars()) 613 614 if nksmall is not None: 615 # if nksmall add bandstructure and dos calculations as well 616 bands_ksampling = aobj.KSampling.path_from_structure(ndivsm=nksmall, structure=structure) 617 dos_ksampling = aobj.KSampling.automatic_density(structure=structure, kppa=2000) 618 nscf_multi[0].set_vars(bands_ksampling.to_abivars()) 619 nscf_multi[0].set_vars({'chksymbreak': 0}) 620 nscf_multi[1].set_vars(dos_ksampling.to_abivars()) 621 nscf_multi[1].set_vars({'chksymbreak': 0}) 622 623 nscf_inputs = nscf_multi.split_datasets() 624 625 # create screening and sigma inputs 626 627 #if scr_nband is None: 628 # scr_nband = nscf_nband_nscf 629 #if sigma_nband is None: 630 # sigma_nband = nscf_nband_nscf 631 632 if 'cd' in response_models: 633 hilbert = aobj.HilbertTransform(nomegasf=100, domegasf=None, spmeth=1, nfreqre=None, freqremax=None, nfreqim=None, 634 freqremin=None) 635 scr_inputs = [] 636 sigma_inputs = [] 637 #print("ecuteps", ecuteps, "nscf_nband", nscf_nband) 638 639 for response_model in response_models: 640 for ecuteps_v in ecuteps: 641 for nscf_nband_v in nscf_nband: 642 scr_nband = nscf_nband_v 643 sigma_nband = nscf_nband_v 644 multi = MultiDataset(structure, pseudos, ndtset=2) 645 multi.set_vars(nscf_ksampling.to_abivars()) 646 multi.set_vars(nscf_electrons.to_abivars()) 647 multi.set_vars(extra_abivars_all) 648 multi.set_vars(extra_abivars_gw) 649 if response_model == 'cd': 650 screening = aobj.Screening(ecuteps_v, scr_nband, w_type="RPA", sc_mode="one_shot", hilbert=hilbert, 651 ecutwfn=None, inclvkb=inclvkb) 652 self_energy = aobj.SelfEnergy("gw", "one_shot", sigma_nband, ecutsigx, screening) 653 else: 654 ppmodel = response_model 655 screening = aobj.Screening(ecuteps_v, scr_nband, w_type="RPA", sc_mode="one_shot", 656 hilbert=None, ecutwfn=None, inclvkb=inclvkb) 657 self_energy = aobj.SelfEnergy("gw", "one_shot", sigma_nband, ecutsigx, screening, 658 gw_qprange=gw_qprange, ppmodel=ppmodel) 659 multi[0].set_vars(screening.to_abivars()) 660 multi[0].set_vars(_stopping_criterion("screening", accuracy)) # Dummy 661 multi[1].set_vars(self_energy.to_abivars()) 662 multi[1].set_vars(_stopping_criterion("sigma", accuracy)) # Dummy 663 664 scr_input, sigma_input = multi.split_datasets() 665 scr_inputs.append(scr_input) 666 sigma_inputs.append(sigma_input) 667 668 return scf_inputs, nscf_inputs, scr_inputs, sigma_inputs 669 670 671def bse_with_mdf_inputs(structure, pseudos, 672 scf_kppa, nscf_nband, nscf_ngkpt, nscf_shiftk, 673 ecuteps, bs_loband, bs_nband, mbpt_sciss, mdf_epsinf, 674 ecut=None, pawecutdg=None, 675 exc_type="TDA", bs_algo="haydock", accuracy="normal", spin_mode="polarized", 676 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None): 677 """ 678 Returns a |MultiDataset| object that performs a GS + NSCF + Bethe-Salpeter calculation. 679 The self-energy corrections are approximated with the scissors operator. 680 The screening is modeled with the model dielectric function. 681 682 Args: 683 structure: |Structure| object. 684 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 685 scf_kppa: Defines the sampling used for the SCF run. 686 nscf_nband: Number of bands included in the NSCF run. 687 nscf_ngkpt: Divisions of the k-mesh used for the NSCF and the BSE run. 688 nscf_shiftk: Shifts used for the NSCF and the BSE run. 689 ecuteps: Cutoff energy [Ha] for the screening matrix. 690 bs_loband: Index of the first occupied band included the e-h basis set 691 (ABINIT convention i.e. first band starts at 1). 692 Can be scalar or array of shape (nsppol,) 693 bs_nband: Highest band idex used for the construction of the e-h basis set. 694 mbpt_sciss: Scissor energy in Hartree. 695 mdf_epsinf: Value of the macroscopic dielectric function used in expression for the model dielectric function. 696 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 697 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized from the pseudos 698 according to accuracy) 699 exc_type: Approximation used for the BSE Hamiltonian (Tamm-Dancoff or coupling). 700 bs_algo: Algorith for the computatio of the macroscopic dielectric function. 701 accuracy: Accuracy of the calculation. 702 spin_mode: Spin polarization. 703 smearing: Smearing technique. 704 charge: Electronic charge added to the unit cell. 705 scf_algorithm: Algorithm used for solving the SCF cycle. 706 """ 707 structure = Structure.as_structure(structure) 708 multi = MultiDataset(structure, pseudos, ndtset=3) 709 710 # Set the cutoff energies. 711 d = _find_ecut_pawecutdg(ecut, pawecutdg, multi.pseudos, accuracy) 712 multi.set_vars(ecut=d.ecut, ecutwfn=d.ecut, pawecutdg=d.pawecutdg) 713 714 # Ground-state 715 scf_ksampling = aobj.KSampling.automatic_density(structure, scf_kppa, chksymbreak=0) 716 717 scf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 718 charge=charge, nband=None, fband=None) 719 720 if scf_electrons.nband is None: 721 scf_electrons.nband = _find_scf_nband(structure, multi.pseudos, scf_electrons) 722 723 multi[0].set_vars(scf_ksampling.to_abivars()) 724 multi[0].set_vars(scf_electrons.to_abivars()) 725 multi[0].set_vars(_stopping_criterion("scf", accuracy)) 726 727 # NSCF calculation with the randomly-shifted k-mesh. 728 nscf_ksampling = aobj.KSampling.monkhorst(nscf_ngkpt, shiftk=nscf_shiftk, chksymbreak=0) 729 730 nscf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm={"iscf": -2}, 731 charge=charge, nband=nscf_nband, fband=None) 732 733 multi[1].set_vars(nscf_ksampling.to_abivars()) 734 multi[1].set_vars(nscf_electrons.to_abivars()) 735 multi[1].set_vars(_stopping_criterion("nscf", accuracy)) 736 737 # BSE calculation. 738 exc_ham = aobj.ExcHamiltonian(bs_loband, bs_nband, mbpt_sciss, coulomb_mode="model_df", ecuteps=ecuteps, 739 spin_mode=spin_mode, mdf_epsinf=mdf_epsinf, exc_type=exc_type, algo=bs_algo, 740 bs_freq_mesh=None, with_lf=True, zcut=None) 741 742 multi[2].set_vars(nscf_ksampling.to_abivars()) 743 multi[2].set_vars(nscf_electrons.to_abivars()) 744 multi[2].set_vars(exc_ham.to_abivars()) 745 #multi[2].set_vars(_stopping_criterion("nscf", accuracy)) 746 747 # TODO: Cannot use istwfk != 1. 748 multi.set_vars(istwfk="*1") 749 750 return multi 751 752 753def scf_phonons_inputs(structure, pseudos, kppa, 754 ecut=None, pawecutdg=None, scf_nband=None, accuracy="normal", spin_mode="polarized", 755 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None): 756 # TODO: Please check the unused variables in the function 757 """ 758 Returns a list of input files for performing phonon calculations. 759 GS input + the input files for the phonon calculation. 760 761 Args: 762 structure: |Structure| object. 763 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 764 kppa: Defines the sampling used for the SCF run. 765 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 766 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized from the 767 pseudos according to accuracy) 768 scf_nband: Number of bands for SCF run. If scf_nband is None, nband is automatically initialized from the list of 769 pseudos, the structure and the smearing option. 770 accuracy: Accuracy of the calculation. 771 spin_mode: Spin polarization. 772 smearing: Smearing technique. 773 charge: Electronic charge added to the unit cell. 774 scf_algorithm: Algorithm used for solving of the SCF cycle. 775 """ 776 # Build the input file for the GS run. 777 gs_inp = AbinitInput(structure=structure, pseudos=pseudos) 778 779 # Set the cutoff energies. 780 gs_inp.set_vars(_find_ecut_pawecutdg(ecut, pawecutdg, gs_inp.pseudos, accuracy)) 781 782 ksampling = aobj.KSampling.automatic_density(gs_inp.structure, kppa, chksymbreak=0) 783 gs_inp.set_vars(ksampling.to_abivars()) 784 gs_inp.set_vars(tolvrs=1.0e-18) 785 786 # Get the qpoints in the IBZ. Note that here we use a q-mesh with ngkpt=(4,4,4) and shiftk=(0,0,0) 787 # i.e. the same parameters used for the k-mesh in gs_inp. 788 qpoints = gs_inp.abiget_ibz(ngkpt=(4, 4, 4), shiftk=(0, 0, 0), kptopt=1).points 789 #print("get_ibz qpoints:", qpoints) 790 791 # Build the input files for the q-points in the IBZ. 792 #ph_inputs = MultiDataset(gs_inp.structure, pseudos=gs_inp.pseudos, ndtset=len(qpoints)) 793 794 ph_inputs = MultiDataset.replicate_input(gs_inp, ndtset=len(qpoints)) 795 796 for ph_inp, qpt in zip(ph_inputs, qpoints): 797 # Response-function calculation for phonons. 798 ph_inp.set_vars( 799 rfphon=1, # Will consider phonon-type perturbation 800 nqpt=1, # One wavevector is to be considered 801 qpt=qpt, # This wavevector is q=0 (Gamma) 802 tolwfr=1.0e-20, 803 kptopt=3, # TODO: One could use symmetries for Gamma. 804 ) 805 #rfatpol 1 1 # Only the first atom is displaced 806 #rfdir 1 0 0 # Along the first reduced coordinate axis 807 #kptopt 2 # Automatic generation of k points, taking 808 809 irred_perts = ph_inp.abiget_irred_phperts() 810 # TODO irred_perts is not used ?? 811 812 #for pert in irred_perts: 813 # #print(pert) 814 # # TODO this will work for phonons, but not for the other types of perturbations. 815 # ph_inp = q_inp.deepcopy() 816 # rfdir = 3 * [0] 817 # rfdir[pert.idir -1] = 1 818 # ph_inp.set_vars( 819 # rfdir=rfdir, 820 # rfatpol=[pert.ipert, pert.ipert] 821 # ) 822 # ph_inputs.append(ph_inp) 823 824 # Split input into gs_inp and ph_inputs 825 all_inps = [gs_inp] 826 all_inps.extend(ph_inputs.split_datasets()) 827 828 return all_inps 829 830 831def phonons_from_gsinput(gs_inp, ph_ngqpt=None, qpoints=None, with_ddk=True, with_dde=True, with_bec=False, 832 ph_tol=None, ddk_tol=None, dde_tol=None, wfq_tol=None, qpoints_to_skip=None, manager=None): 833 """ 834 Returns a list of inputs in the form of a MultiDataset to perform phonon calculations, based on 835 a ground state |AbinitInput|. 836 It will determine if WFQ files should be calculated for some q points and add the NSCF AbinitInputs to the set. 837 The inputs have the following tags, according to their function: "ddk", "dde", "nscf", "ph_q_pert". 838 All of them have the tag "phonon". 839 840 Args: 841 gs_inp: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK. 842 ph_ngqpt: a list of three integers representing the gamma centered q-point grid used for the calculation. 843 If None and qpoint==None the ngkpt value present in the gs_input will be used. 844 Incompatible with qpoints. 845 qpoints: a list of coordinates of q points in reduced coordinates for which the phonon perturbations will 846 be calculated. Incompatible with ph_ngqpt. 847 with_ddk: If True, if Gamma is included in the list of qpoints it will add inputs for the calculations of 848 the DDK. 849 with_dde: If True, if Gamma is included in the list of qpoints it will add inputs for the calculations of 850 the DDE. Automatically sets with_ddk=True. 851 with_bec: If Truem if Gamma is included in the list of qpoints the DDE will be calculated in the same 852 input as the phonons. This will allow to determine the BECs. 853 Automatically sets with_ddk=True and with_dde=False. 854 ph_tol: a dictionary with a single key defining the type of tolerance used for the phonon calculations and 855 its value. Default: {"tolvrs": 1.0e-10}. 856 ddk_tol: a dictionary with a single key defining the type of tolerance used for the DDK calculations and 857 its value. Default: {"tolwfr": 1.0e-22}. 858 dde_tol: a dictionary with a single key defining the type of tolerance used for the DDE calculations and 859 its value. Default: {"tolvrs": 1.0e-10}. 860 wfq_tol: a dictionary with a single key defining the type of tolerance used for the NSCF calculations of 861 the WFQ and its value. Default {"tolwfr": 1.0e-22}. 862 qpoints_to_skip: a list of coordinates of q points in reduced coordinates that will be skipped. 863 Useful when calculating multiple grids for the same system to avoid duplicate calculations. 864 If a DDB needs to be extended with more q points use e.g. ddb.qpoints.to_array(). 865 manager: |TaskManager| of the task. If None, the manager is initialized from the config file. 866 """ 867 gs_inp = gs_inp.deepcopy() 868 gs_inp.pop_irdvars() 869 870 if with_dde: 871 with_ddk = True 872 873 if with_bec: 874 with_ddk = True 875 with_dde = False 876 877 if ph_tol is None: 878 ph_tol = {"tolvrs": 1.0e-10} 879 880 if ddk_tol is None: 881 ddk_tol = {"tolwfr": 1.0e-22} 882 883 if dde_tol is None: 884 dde_tol = {"tolvrs": 1.0e-10} 885 886 if wfq_tol is None: 887 wfq_tol = {"tolwfr": 1.0e-22} 888 889 multi = [] 890 891 if qpoints is not None and ph_ngqpt is not None: 892 raise ValueError("ph_ngqpt and qpoints can't be used together") 893 894 if qpoints is None: 895 if ph_ngqpt is None: 896 ph_ngqpt = np.array(gs_inp["ngkpt"]) 897 else: 898 ph_ngqpt = np.array(ph_ngqpt) 899 900 qpoints = gs_inp.abiget_ibz(ngkpt=ph_ngqpt, shiftk=(0, 0, 0), kptopt=1, manager=manager).points 901 902 if qpoints_to_skip: 903 preserved_qpoints = [] 904 for q in qpoints: 905 if not any(np.allclose(q, ddb_q) for ddb_q in qpoints_to_skip): 906 preserved_qpoints.append(q) 907 qpoints = np.array(preserved_qpoints) 908 909 if ph_ngqpt is None or any(gs_inp["ngkpt"] % ph_ngqpt != 0): 910 # find which q points are needed and build nscf inputs to calculate the WFQ 911 kpts = gs_inp.abiget_ibz(shiftk=(0, 0, 0), kptopt=3, manager=manager).points.tolist() 912 nscf_qpt = [] 913 for q in qpoints: 914 if list(q) not in kpts: 915 nscf_qpt.append(q) 916 if nscf_qpt: 917 multi_nscf = MultiDataset.replicate_input(gs_inp, len(nscf_qpt)) 918 multi_nscf.set_vars(kptopt=3, nqpt=1, iscf=-2) 919 if wfq_tol: 920 multi_nscf.set_vars(**wfq_tol) 921 else: 922 multi_nscf.set_vars(tolwfr=1e-22) 923 for q, nscf_inp in zip(nscf_qpt, multi_nscf): 924 nscf_inp.set_vars(qpt=q) 925 926 multi_nscf.add_tags(atags.NSCF) 927 928 multi.extend(multi_nscf) 929 930 # Build the input files for the q-points in the IBZ. 931 # Response-function calculation for phonons. 932 for qpt in qpoints: 933 if np.allclose(qpt, 0): 934 if with_ddk: 935 multi_ddk = gs_inp.make_ddk_inputs(tolerance=ddk_tol) 936 multi_ddk.add_tags(atags.DDK) 937 multi.extend(multi_ddk) 938 if with_dde: 939 multi_dde = gs_inp.make_dde_inputs(dde_tol, manager=manager) 940 multi_dde.add_tags(atags.DDE) 941 multi.extend(multi_dde) 942 elif with_bec: 943 multi_bec = gs_inp.make_bec_inputs(ph_tol, manager=manager) 944 multi_bec.add_tags(atags.BEC) 945 multi.extend(multi_bec) 946 continue 947 948 multi_ph_q = gs_inp.make_ph_inputs_qpoint(qpt, ph_tol) 949 multi_ph_q.add_tags(atags.PH_Q_PERT) 950 multi.extend(multi_ph_q) 951 952 multi = MultiDataset.from_inputs(multi) 953 multi.add_tags(atags.PHONON) 954 955 return multi 956 957 958def piezo_elastic_inputs_from_gsinput(gs_inp, ddk_tol=None, rf_tol=None, ddk_split=False, rf_split=False, 959 manager=None): 960 """ 961 Returns a |MultiDataset| for performing elastic and piezoelectric constants calculations. 962 GS input + the input files for the elastic and piezoelectric constants calculation. 963 964 Args: 965 gs_inp: Ground State input to build piezo elastic inputs from. 966 ddk_tol: Tolerance for the DDK calculation (i.e. {"tolwfr": 1.0e-20}). 967 rf_tol: Tolerance for the Strain RF calculations (i.e. {"tolvrs": 1.0e-12}). 968 ddk_split: Whether to split the DDK calculations. 969 rf_split: whether to split the RF calculations. 970 manager: |TaskManager| of the task. If None, the manager is initialized from the config file. 971 """ 972 # Ddk input(s) 973 if ddk_split: 974 multi = gs_inp.make_ddk_inputs(tolerance=ddk_tol) 975 else: 976 ddk_inp = gs_inp.deepcopy() 977 978 ddk_inp.set_vars( 979 rfelfd=2, # Activate the calculation of the d/dk perturbation 980 rfdir=(1,1,1), # All directions 981 nqpt=1, # One wavevector is to be considered 982 qpt=(0, 0, 0), # q-wavevector. 983 kptopt=3, # Take into account time-reversal symmetry. 984 iscf=-3, # The d/dk perturbation must be treated in a non-self-consistent way 985 paral_kgb=0 986 ) 987 if ddk_tol is None: 988 ddk_tol = {"tolwfr": 1.0e-20} 989 990 if len(ddk_tol) != 1 or any(k not in _tolerances for k in ddk_tol): 991 raise ValueError("Invalid tolerance: {}".format(ddk_tol)) 992 ddk_inp.pop_tolerances() 993 ddk_inp.set_vars(ddk_tol) 994 # Adding buffer to help convergence ... 995 if 'nbdbuf' not in ddk_inp: 996 nbdbuf = max(int(0.1*ddk_inp['nband']), 4) 997 ddk_inp.set_vars(nband=ddk_inp['nband']+nbdbuf, nbdbuf=nbdbuf) 998 999 multi = MultiDataset.from_inputs([ddk_inp]) 1000 multi.add_tags(atags.DDK) 1001 1002 # Response Function input(s) 1003 if rf_split: 1004 multi_rf = gs_inp.make_strain_perts_inputs(tolerance=rf_tol, manager=manager) 1005 else: 1006 rf_inp = gs_inp.deepcopy() 1007 1008 rf_inp.set_vars(rfphon=1, # Atomic displacement perturbation 1009 rfatpol=(1,len(gs_inp.structure)), # Perturbation of all atoms 1010 rfstrs=3, # Do the strain perturbations 1011 rfdir=(1,1,1), # All directions 1012 nqpt=1, # One wavevector is to be considered 1013 qpt=(0, 0, 0), # q-wavevector. 1014 kptopt=3, # Take into account time-reversal symmetry. 1015 iscf=7, # The rfstrs perturbation must be treated in a 1016 # self-consistent way 1017 paral_kgb=0 1018 ) 1019 1020 if rf_tol is None: 1021 rf_tol = {"tolvrs": 1.0e-12} 1022 1023 if len(rf_tol) != 1 or any(k not in _tolerances for k in rf_tol): 1024 raise ValueError("Invalid tolerance: {}".format(rf_tol)) 1025 rf_inp.pop_tolerances() 1026 rf_inp.set_vars(rf_tol) 1027 1028 # Adding buffer to help convergence ... 1029 if 'nbdbuf' not in rf_inp: 1030 nbdbuf = max(int(0.1*rf_inp['nband']), 4) 1031 rf_inp.set_vars(nband=rf_inp['nband']+nbdbuf, nbdbuf=nbdbuf) 1032 1033 multi_rf = MultiDataset.from_inputs([rf_inp]) 1034 multi_rf.add_tags([atags.DFPT, atags.STRAIN]) 1035 for inp in multi_rf: 1036 if inp.get('rfphon', 0) == 1: 1037 inp.add_tags(atags.PHONON) 1038 1039 multi.extend(multi_rf) 1040 1041 return multi 1042 1043 1044def scf_piezo_elastic_inputs(structure, pseudos, kppa, ecut=None, pawecutdg=None, scf_nband=None, 1045 accuracy="normal", spin_mode="polarized", 1046 smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None, 1047 ddk_tol=None, rf_tol=None, ddk_split=False, rf_split=False): 1048 1049 """ 1050 Returns a |MultiDataset| for performing elastic and piezoelectric constants calculations. 1051 GS input + the input files for the elastic and piezoelectric constants calculation. 1052 1053 Args: 1054 structure: |Structure| object. 1055 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 1056 kppa: Defines the sampling used for the SCF run. 1057 ecut: cutoff energy in Ha (if None, ecut is initialized from the pseudos according to accuracy) 1058 pawecutdg: cutoff energy in Ha for PAW double-grid (if None, pawecutdg is initialized from the 1059 pseudos according to accuracy) 1060 scf_nband: Number of bands for SCF run. If scf_nband is None, nband is automatically initialized 1061 from the list of pseudos, the structure and the smearing option. 1062 accuracy: Accuracy of the calculation. 1063 spin_mode: Spin polarization. 1064 smearing: Smearing technique. 1065 charge: Electronic charge added to the unit cell. 1066 scf_algorithm: Algorithm used for solving of the SCF cycle. 1067 ddk_tol: Tolerance for the Ddk calculation (i.e. {"tolwfr": 1.0e-20}). 1068 rf_tol: Tolerance for the Strain RF calculations (i.e. {"tolvrs": 1.0e-12}). 1069 ddk_split: Whether to split the ddk calculations. 1070 rf_split: whether to split the RF calculations. 1071 """ 1072 # Build the input file for the GS run. 1073 gs_inp = scf_input(structure=structure, pseudos=pseudos, kppa=kppa, ecut=ecut, pawecutdg=pawecutdg, 1074 nband=scf_nband, accuracy=accuracy, spin_mode=spin_mode, smearing=smearing, charge=charge, 1075 scf_algorithm=scf_algorithm, shift_mode="Gamma-centered") 1076 1077 # Adding buffer to help convergence ... 1078 nbdbuf = max(int(0.1*gs_inp['nband']), 4) 1079 gs_inp.set_vars(nband=gs_inp['nband']+nbdbuf, nbdbuf=nbdbuf) 1080 1081 multi = MultiDataset.from_inputs([gs_inp]) 1082 1083 piezo_elastic_inputs = piezo_elastic_inputs_from_gsinput(gs_inp=gs_inp, ddk_tol=ddk_tol, rf_tol=rf_tol) 1084 1085 multi.extend(piezo_elastic_inputs) 1086 1087 return multi 1088 1089 1090def scf_input(structure, pseudos, kppa=None, ecut=None, pawecutdg=None, nband=None, accuracy="normal", 1091 spin_mode="polarized", smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None, 1092 shift_mode="Monkhorst-Pack"): 1093 """ 1094 Returns an |AbinitInput| object for standard GS calculations. 1095 """ 1096 structure = Structure.as_structure(structure) 1097 1098 abinit_input = AbinitInput(structure, pseudos) 1099 1100 # Set the cutoff energies. 1101 abinit_input.set_vars(_find_ecut_pawecutdg(ecut, pawecutdg, abinit_input.pseudos, accuracy)) 1102 1103 # SCF calculation. 1104 kppa = _DEFAULTS.get("kppa") if kppa is None else kppa 1105 shift_mode = ShiftMode.from_object(shift_mode) 1106 shifts = _get_shifts(shift_mode, structure) 1107 scf_ksampling = aobj.KSampling.automatic_density(structure, kppa, chksymbreak=0, shifts=shifts) 1108 scf_electrons = aobj.Electrons(spin_mode=spin_mode, smearing=smearing, algorithm=scf_algorithm, 1109 charge=charge, nband=nband, fband=None) 1110 1111 if spin_mode == "polarized": 1112 abinit_input.set_autospinat() 1113 1114 if scf_electrons.nband is None: 1115 scf_electrons.nband = _find_scf_nband(structure, abinit_input.pseudos, scf_electrons, 1116 abinit_input.get('spinat', None)) 1117 1118 abinit_input.set_vars(scf_ksampling.to_abivars()) 1119 abinit_input.set_vars(scf_electrons.to_abivars()) 1120 abinit_input.set_vars(_stopping_criterion("scf", accuracy)) 1121 1122 return abinit_input 1123 1124 1125def ebands_from_gsinput(gsinput, nband=None, ndivsm=15, accuracy="normal"): 1126 """ 1127 Return an |AbinitInput| object to compute a band structure from a GS SCF input. 1128 1129 Args: 1130 gsinput: 1131 nband: 1132 ndivsm: 1133 accuracy: 1134 1135 Return: |AbinitInput| 1136 """ 1137 # create a copy to avoid messing with the previous input 1138 bands_input = gsinput.deepcopy() 1139 1140 bands_input.pop_irdvars() 1141 1142 nscf_ksampling = aobj.KSampling.path_from_structure(ndivsm, gsinput.structure) 1143 if nband is None: 1144 nband = gsinput.get("nband", gsinput.structure.num_valence_electrons(gsinput.pseudos)) + 10 1145 1146 bands_input.set_vars(nscf_ksampling.to_abivars()) 1147 bands_input.set_vars(nband=nband, iscf=-2) 1148 bands_input.set_vars(_stopping_criterion("nscf", accuracy)) 1149 1150 return bands_input 1151 1152 1153def dos_from_gsinput(gsinput, dos_kppa, nband=None, accuracy="normal", pdos=False): 1154 1155 # create a copy to avoid messing with the previous input 1156 dos_input = gsinput.deepcopy() 1157 dos_input.pop_irdvars() 1158 1159 dos_ksampling = aobj.KSampling.automatic_density(dos_input.structure, dos_kppa, chksymbreak=0) 1160 dos_input.set_vars(dos_ksampling.to_abivars()) 1161 dos_input.set_vars(iscf=-2, ionmov=0) 1162 dos_input.set_vars(_stopping_criterion("nscf", accuracy)) 1163 1164 if pdos: 1165 # FIXME 1166 raise NotImplementedError() 1167 1168 return dos_input 1169 1170 1171def ioncell_relax_from_gsinput(gsinput, accuracy="normal"): 1172 1173 ioncell_input = gsinput.deepcopy() 1174 ioncell_input.pop_irdvars() 1175 1176 ioncell_relax = aobj.RelaxationMethod.atoms_and_cell(atoms_constraints=None) 1177 ioncell_input.set_vars(ioncell_relax.to_abivars()) 1178 ioncell_input.set_vars(_stopping_criterion("relax", accuracy)) 1179 1180 return ioncell_input 1181 1182 1183def hybrid_oneshot_input(gsinput, functional="hse06", ecutsigx=None, gw_qprange=1): 1184 1185 hybrid_input = gsinput.deepcopy() 1186 hybrid_input.pop_irdvars() 1187 1188 functional = functional.lower() 1189 if functional == 'hse06': 1190 gwcalctyp = 115 1191 icutcoul = 5 1192 rcut = 9.090909 1193 elif functional == 'pbe0': 1194 gwcalctyp = 215 1195 icutcoul = 6 1196 rcut = 0. 1197 elif functional == 'b3lyp': 1198 gwcalctyp = 315 1199 icutcoul = 6 1200 rcut = 0. 1201 else: 1202 raise ValueError("Unknow functional {0}.".format(functional)) 1203 1204 ecut = hybrid_input['ecut'] 1205 ecutsigx = ecutsigx or 2*ecut 1206 1207 hybrid_input.set_vars(optdriver=4, gwcalctyp=gwcalctyp, gw_nstep=1, gwpara=2, icutcoul=icutcoul, rcut=rcut, 1208 gw_qprange=gw_qprange, ecutwfn=ecut*0.995, ecutsigx=ecutsigx) 1209 1210 return hybrid_input 1211 1212 1213def hybrid_scf_input(gsinput, functional="hse06", ecutsigx=None, gw_qprange=1): 1214 1215 hybrid_input = hybrid_oneshot_input(gsinput=gsinput, functional=functional, ecutsigx=ecutsigx, gw_qprange=gw_qprange) 1216 1217 hybrid_input['gwcalctyp'] += 10 1218 1219 return hybrid_input 1220 1221 1222def scf_for_phonons(structure, pseudos, kppa=None, ecut=None, pawecutdg=None, nband=None, accuracy="normal", 1223 spin_mode="polarized", smearing="fermi_dirac:0.1 eV", charge=0.0, scf_algorithm=None, 1224 shift_mode="Symmetric"): 1225 1226 abiinput = scf_input(structure=structure, pseudos=pseudos, kppa=kppa, ecut=ecut, pawecutdg=pawecutdg, nband=nband, 1227 accuracy=accuracy, spin_mode=spin_mode, smearing=smearing, charge=charge, 1228 scf_algorithm=scf_algorithm, shift_mode=shift_mode) 1229 1230 nbdbuf = 4 1231 # with no smearing set the minimum number of bands plus some nbdbuf 1232 if smearing is None: 1233 nval = structure.num_valence_electrons(pseudos) 1234 nval -= abiinput['charge'] 1235 nband = int(round(nval / 2) + nbdbuf) 1236 abiinput.set_vars(nband=nband) 1237 1238 # enforce symmetries and add a buffer of bands to ease convergence with tolwfr 1239 abiinput.set_vars(chksymbreak=1, nbdbuf=nbdbuf, tolwfr=1.e-22) 1240 1241 return abiinput 1242 1243 1244def dte_from_gsinput(gs_inp, use_phonons=True, ph_tol=None, ddk_tol=None, dde_tol=None, 1245 skip_dte_permutations=False, manager=None): 1246 """ 1247 Returns a list of inputs in the form of a |MultiDataset| to perform calculations of non-linear properties, based on 1248 a ground state AbinitInput. 1249 1250 The inputs have the following tags, according to their function: "ddk", "dde", "ph_q_pert" and "dte". 1251 All of them have the tag "dfpt". 1252 1253 Args: 1254 gs_inp: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK. 1255 use_phonons: determine wether the phonon perturbations at gamma should be included or not 1256 ph_tol: a dictionary with a single key defining the type of tolerance used for the phonon calculations and 1257 its value. Default: {"tolvrs": 1.0e-22}. 1258 ddk_tol: a dictionary with a single key defining the type of tolerance used for the DDK calculations and 1259 its value. Default: {"tolwfr": 1.0e-22}. 1260 dde_tol: a dictionary with a single key defining the type of tolerance used for the DDE calculations and 1261 its value. Default: {"tolvrs": 1.0e-22}. 1262 skip_dte_permutations: Since the current version of abinit always performs all the permutations of the 1263 perturbations, even if only one is asked, if True avoids the creation of inputs that will produce 1264 duplicated outputs. 1265 manager: |TaskManager| of the task. If None, the manager is initialized from the config file. 1266 """ 1267 gs_inp = gs_inp.deepcopy() 1268 gs_inp.pop_irdvars() 1269 1270 if ph_tol is None: 1271 ph_tol = {"tolvrs": 1.0e-22} 1272 1273 if ddk_tol is None: 1274 ddk_tol = {"tolwfr": 1.0e-22} 1275 1276 if dde_tol is None: 1277 dde_tol = {"tolvrs": 1.0e-22} 1278 1279 multi = [] 1280 1281 multi_ddk = gs_inp.make_ddk_inputs(tolerance=ddk_tol) 1282 multi_ddk.add_tags(atags.DDK) 1283 multi.extend(multi_ddk) 1284 multi_dde = gs_inp.make_dde_inputs(dde_tol, use_symmetries=False, manager=manager) 1285 multi_dde.add_tags(atags.DDE) 1286 multi.extend(multi_dde) 1287 1288 if use_phonons: 1289 multi_ph = gs_inp.make_ph_inputs_qpoint([0,0,0], ph_tol, manager=manager) 1290 multi_ph.add_tags(atags.PH_Q_PERT) 1291 multi.extend(multi_ph) 1292 1293 # non-linear calculations do not accept more bands than those in the valence. Set the correct values. 1294 # Do this as last, so not to interfere with the the generation of the other steps. 1295 nval = gs_inp.structure.num_valence_electrons(gs_inp.pseudos) 1296 nval -= gs_inp['charge'] 1297 nband = int(round(nval / 2)) 1298 gs_inp.set_vars(nband=nband) 1299 gs_inp.pop('nbdbuf', None) 1300 multi_dte = gs_inp.make_dte_inputs(phonon_pert=use_phonons, skip_permutations=skip_dte_permutations, 1301 manager=manager) 1302 multi_dte.add_tags(atags.DTE) 1303 multi.extend(multi_dte) 1304 1305 multi = MultiDataset.from_inputs(multi) 1306 multi.add_tags(atags.DFPT) 1307 1308 return multi 1309 1310 1311def dfpt_from_gsinput(gs_inp, ph_ngqpt=None, qpoints=None, do_ddk=True, do_dde=True, do_strain=True, 1312 do_dte=False, ph_tol=None, ddk_tol=None, dde_tol=None, wfq_tol=None, strain_tol=None, 1313 skip_dte_permutations=False, manager=None): 1314 """ 1315 Returns a list of inputs in the form of a MultiDataset to perform a set of calculations based on DFPT including 1316 phonons, elastic and non-linear properties. Requires a ground state |AbinitInput| as a starting point. 1317 1318 It will determine if WFQ files should be calculated for some q points and add the NSCF AbinitInputs to the set. 1319 The original input is included and the inputs have the following tags, according to their function: 1320 "scf", "ddk", "dde", "nscf", "ph_q_pert", "strain", "dte", "dfpt". 1321 1322 N.B. Currently (version 8.8.3) anaddb does not support a DDB containing both 2nd order derivatives with qpoints 1323 different from gamma AND 3rd oreder derivatives. The calculations could be run, but the global DDB will not 1324 be directly usable as is. 1325 1326 Args: 1327 gs_inp: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK. 1328 ph_ngqpt: a list of three integers representing the gamma centered q-point grid used for the calculation. 1329 If None and qpoint==None the ngkpt value present in the gs_input will be used. 1330 Incompatible with qpoints. 1331 qpoints: a list of coordinates of q points in reduced coordinates for which the phonon perturbations will 1332 be calculated. Incompatible with ph_ngqpt. 1333 do_ddk: If True, if Gamma is included in the list of qpoints it will add inputs for the calculations of 1334 the DDK. 1335 do_dde: If True, if Gamma is included in the list of qpoints it will add inputs for the calculations of 1336 the DDE. Automatically sets with_ddk=True. 1337 do_strain: If True inputs for the strain perturbations will be included. 1338 do_dte: If True inputs for the non-linear perturbations will be included. The phonon non-linear perturbations 1339 will be included only if a phonon calculation at gamma is present. The caller is responsible for 1340 adding it. Automatically sets with_dde=True. 1341 ph_tol: a dictionary with a single key defining the type of tolerance used for the phonon calculations and 1342 its value. Default: {"tolvrs": 1.0e-10}. 1343 ddk_tol: a dictionary with a single key defining the type of tolerance used for the DDK calculations and 1344 its value. Default: {"tolwfr": 1.0e-22}. 1345 dde_tol: a dictionary with a single key defining the type of tolerance used for the DDE calculations and 1346 its value. Default: {"tolvrs": 1.0e-10}. 1347 wfq_tol: a dictionary with a single key defining the type of tolerance used for the NSCF calculations of 1348 the WFQ and its value. Default {"tolwfr": 1.0e-22}. 1349 strain_tol: dictionary with a single key defining the type of tolerance used for the strain calculations of 1350 and its value. Default {"tolvrs": 1.0e-12}. 1351 skip_dte_permutations: Since the current version of abinit always performs all the permutations of the 1352 perturbations, even if only one is asked, if True avoids the creation of inputs that will produce 1353 duplicated outputs. 1354 manager: |TaskManager| of the task. If None, the manager is initialized from the config file. 1355 """ 1356 1357 if ph_tol is None: 1358 ph_tol = {"tolvrs": 1.0e-10} 1359 if ddk_tol is None: 1360 ddk_tol = {"tolwfr": 1.0e-22} 1361 if dde_tol is None: 1362 dde_tol = {"tolvrs": 1.0e-10} 1363 if wfq_tol is None: 1364 wfq_tol = {"tolwfr": 1.0e-22} 1365 if strain_tol is None: 1366 strain_tol = {"tolvrs": 1.0e-12} 1367 1368 if do_dde: 1369 do_ddk = True 1370 1371 if do_dte: 1372 do_dde = True 1373 1374 multi = MultiDataset.from_inputs([gs_inp]) 1375 multi[0].add_tags(atags.SCF) 1376 1377 do_phonons = ph_ngqpt is not None or qpoints is not None 1378 has_gamma = False 1379 if do_phonons: 1380 multi.extend(phonons_from_gsinput(gs_inp, ph_ngqpt=ph_ngqpt, qpoints=qpoints, with_ddk=False, with_dde=False, 1381 with_bec=False, ph_tol=ph_tol, ddk_tol=ddk_tol, dde_tol=dde_tol, 1382 wfq_tol=wfq_tol, qpoints_to_skip=None, manager=manager)) 1383 has_gamma = ph_ngqpt is not None or any(np.allclose(q, [0, 0, 0]) for q in qpoints) 1384 1385 if do_ddk: 1386 multi_ddk = gs_inp.make_ddk_inputs(tolerance=ddk_tol) 1387 multi_ddk.add_tags(atags.DDK) 1388 multi.extend(multi_ddk) 1389 if do_dde: 1390 multi_dde = gs_inp.make_dde_inputs(dde_tol, use_symmetries=not do_dte, manager=manager) 1391 multi_dde.add_tags(atags.DDE) 1392 multi.extend(multi_dde) 1393 1394 if do_strain: 1395 multi_strain = gs_inp.make_strain_perts_inputs(tolerance=strain_tol, manager=manager, phonon_pert=False, 1396 kptopt=2) 1397 multi_strain.add_tags([atags.DFPT, atags.STRAIN]) 1398 multi.extend(multi_strain) 1399 1400 if do_dte: 1401 # non-linear calculations do not accept more bands than those in the valence. Set the correct values. 1402 nval = gs_inp.structure.num_valence_electrons(gs_inp.pseudos) 1403 nval -= gs_inp['charge'] 1404 nband = int(round(nval / 2)) 1405 gs_inp_copy = gs_inp.deepcopy() 1406 gs_inp_copy.set_vars(nband=nband) 1407 gs_inp_copy.pop('nbdbuf', None) 1408 multi_dte = gs_inp_copy.make_dte_inputs(phonon_pert=do_phonons and has_gamma, 1409 skip_permutations=skip_dte_permutations, manager=manager) 1410 multi_dte.add_tags([atags.DTE, atags.DFPT]) 1411 multi.extend(multi_dte) 1412 1413 return multi 1414 1415 1416def conduc_from_inputs(scf_input, nscf_input, tmesh, ddb_ngqpt, eph_ngqpt_fine, sigma_erange, boxcutmin=1.1, mixprec=1): 1417 """ 1418 Returns a list of inputs in the form of a MultiDataset to perform a set of calculations to determine conductivity. 1419 This part require a ground state |AbinitInput| and a non self-consistent |AbinitInput|. You will also need 1420 a work to get DDB and DVDB since |ConducWork| needs these files. 1421 1422 Args: 1423 scf_input: |AbinitInput| representing a ground state calculation, the SCF performed to get the WFK. 1424 nscf_input: |AbinitInput| representing a nscf ground state calculation, the NSCF performed to get the WFK. 1425 most parameters for subsequent tasks will be taken from this inputs. 1426 tmesh: The mesh of temperature (in Kelvin) where we calculate the conductivity. 1427 ddb_ngqpt: the coarse grid of q-points used to compute the DDB and DVDB files in the previous phonon_work. 1428 eph_ngqpt_fine: the fine grid of q-points used for the Fourier nterpolation. 1429 boxcutmin: For the last task only, 1.1 is often used to decrease memory and is faster over the Abinit default of 2. 1430 mixprec: For the last task only, 1 is often used to make the EPH calculation faster. Note that Abinit default is 0. 1431 """ 1432 # Create a MultiDataset from scf input 1433 multi = MultiDataset.from_inputs([scf_input]) 1434 1435 # Add 2 times the nscf input at the end of the MultiDataset 1436 extension = MultiDataset.replicate_input(nscf_input, 2) 1437 multi.extend(extension) 1438 1439 # Modify the second nscf input to get a task that interpolate the DVDB 1440 #multi[2].pop_vars("iscf") 1441 #multi[2].set_vars(irdden=0, optdriver=7, 1442 # ddb_ngqpt=ddb_ngqpt, 1443 # eph_task=5, 1444 # eph_ngqpt_fine=eph_ngqpt_fine) 1445 1446 # Modify the third nscf input to get a conductivity task 1447 multi[2].pop_vars("iscf") 1448 multi[2].set_vars(irdden=0, 1449 optdriver=7, 1450 ddb_ngqpt=ddb_ngqpt, 1451 eph_ngqpt_fine=eph_ngqpt_fine, 1452 eph_task=-4, 1453 tmesh=tmesh, 1454 sigma_erange=sigma_erange, 1455 boxcutmin=boxcutmin, 1456 mixprec=mixprec) 1457 1458 return multi 1459 1460 1461def conduc_kerange_from_inputs(scf_input, nscf_input, tmesh, ddb_ngqpt, eph_ngqpt_fine, 1462 sigma_ngkpt, sigma_erange, einterp=(1, 5, 0, 0), boxcutmin=1.1, mixprec=1): 1463 """ 1464 Returns a list of inputs in the form of a MultiDataset to perform a set of calculations to determine the conductivity. 1465 This part require a ground state |AbinitInput| and a non self-consistent |AbinitInput|. You will also need 1466 a work to get DDB and DVDB since |ConducWork| needs these files. 1467 1468 Args: 1469 scf_input: |AbinitInput| representing a ground state calculation, the SCF performed to get the WFK. 1470 nscf_input: |AbinitInput| representing a nscf ground state calculation, the NSCF performed to get the WFK. 1471 most parameters for subsequent tasks will be taken from this inputs. 1472 ddb_ngqpt: the coarse q-point grid used to get the DDB and DVDB files. 1473 eph_ngqpt_fine: the fine qpoints grid that will be interpolated. 1474 sigma_ngkpt: The fine grid of kpt inside the sigma interval 1475 sigma_erange: The energy range for Sigma_nk 1476 einterp: The interpolation used. By default it is a star-function interpolation. 1477 boxcutmin: For the last task only, 1.1 is often used to decrease memory and is faster over the Abinit default of 2. 1478 mixprec: For the last task only, 1 is often used to make the EPH calculation faster. Note that Abinit default is 0. 1479 """ 1480 # Create a MultiDataset from scf input 1481 multi = MultiDataset.from_inputs([scf_input]) 1482 1483 # Add 4 times the nscf input at the end of the MultiDataset 1484 extension = MultiDataset.replicate_input(nscf_input, 4) 1485 multi.extend(extension) 1486 1487 # Modify the second nscf input to get a task that calculate the kpt in the sigma interval (Kerange.nc file) 1488 multi[2].set_vars(optdriver=8, wfk_task='"wfk_kpts_erange"', kptopt=1, 1489 sigma_ngkpt=sigma_ngkpt, einterp=einterp, sigma_erange=sigma_erange) 1490 1491 # Modify the third nscf input to get a task that add the kpt of Kerange.nc to the WFK file 1492 multi[3].set_vars(optdriver=0, iscf=-2, kptopt=0, ddb_ngqpt=ddb_ngqpt) 1493 1494 # Modify the fourth nscf input to get a task that interpolate the DVDB 1495 #multi[4].pop_vars("iscf") 1496 #multi[4].set_vars(irdden=0, optdriver=7, 1497 # ddb_ngqpt=ddb_ngqpt, 1498 # eph_task=5, 1499 # eph_ngqpt_fine=eph_ngqpt_fine) 1500 1501 # Modify the third nscf input to get a conductivity task 1502 multi[4].pop_vars("iscf") 1503 multi[4].set_vars(irdden=0, 1504 optdriver=7, 1505 ddb_ngqpt=ddb_ngqpt, 1506 eph_ngqpt_fine=eph_ngqpt_fine, 1507 eph_task=-4, 1508 tmesh=tmesh, 1509 sigma_erange=sigma_erange, 1510 ngkpt=sigma_ngkpt, 1511 boxcutmin=boxcutmin, 1512 mixprec=mixprec) 1513 1514 return multi 1515 1516 1517def minimal_scf_input(structure, pseudos): 1518 """ 1519 Provides an input for a calculation with the minimum possible requirements. 1520 Can be used to execute abinit with minimal requirements when needing files 1521 that are produced only after a full calculation completes. 1522 In general this will contain 1 kpt, 1 band, very low cutoff, no polarization, 1523 no smearing. Disables checks on primitive cell and symmetries. 1524 Even for large system it will require small memory allocations and few seconds 1525 to execute. 1526 1527 Args: 1528 structure: |Structure| object. 1529 pseudos: List of filenames or list of |Pseudo| objects or |PseudoTable| object. 1530 """ 1531 1532 inp = scf_input(structure, pseudos, smearing=None, spin_mode="unpolarized") 1533 inp["ngkpt"] = [1, 1, 1] 1534 inp["nshiftk"] = 1 1535 inp["shiftk"] = [[0, 0, 0]] 1536 inp["nstep"] = 0 1537 inp["ecut"] = 3 # should be reasonable, otherwise abinit raises an error 1538 inp["nband"] = 1 1539 inp["chkprim"] = 0 1540 inp["chksymbreak"] = 0 1541 inp["charge"] = structure.num_valence_electrons(inp.pseudos) - 1 1542 inp["boxcutmin"] = 1.2 1543 return inp 1544 1545 1546#FIXME if the pseudos are passed as a PseudoTable the whole table will be serialized, 1547# it would be better to filter on the structure elements 1548class InputFactory(MSONable): 1549 factory_function = None 1550 input_required = True 1551 1552 def __init__(self, *args, **kwargs): 1553 if self.factory_function is None: 1554 raise NotImplementedError('The factory function should be specified') 1555 1556 self.args = args 1557 self.kwargs = kwargs 1558 1559 def build_input(self, previous_input=None): 1560 # make a copy to pop additional parameteres 1561 kwargs = dict(self.kwargs) 1562 decorators = kwargs.pop('decorators', []) 1563 if not isinstance(decorators, (list, tuple)): 1564 decorators = [decorators] 1565 extra_abivars = kwargs.pop('extra_abivars', {}) 1566 if self.input_required: 1567 if not previous_input: 1568 raise ValueError('An input is required for factory function {0}.'.format(self.factory_function.__name__)) 1569 abiinput = self.factory_function(previous_input, *self.args, **kwargs) 1570 else: 1571 abiinput = self.factory_function(*self.args, **kwargs) 1572 1573 for d in decorators: 1574 abiinput = d(abiinput) 1575 abiinput.set_vars(extra_abivars) 1576 1577 return abiinput 1578 1579 @pmg_serialize 1580 def as_dict(self): 1581 # sanitize to avoid numpy arrays and serialize MSONable objects 1582 return jsanitize(dict(args=self.args, kwargs=self.kwargs), strict=True) 1583 1584 @classmethod 1585 def from_dict(cls, d): 1586 dec = MontyDecoder() 1587 return cls(*dec.process_decoded(d['args']), **dec.process_decoded(d['kwargs'])) 1588 1589 1590class BandsFromGsFactory(InputFactory): 1591 factory_function = staticmethod(ebands_from_gsinput) 1592 1593 1594class IoncellRelaxFromGsFactory(InputFactory): 1595 factory_function = staticmethod(ioncell_relax_from_gsinput) 1596 1597 1598class HybridOneShotFromGsFactory(InputFactory): 1599 factory_function = staticmethod(hybrid_oneshot_input) 1600 1601 1602class HybridScfFromGsFactory(InputFactory): 1603 factory_function = staticmethod(hybrid_scf_input) 1604 1605 1606class ScfFactory(InputFactory): 1607 factory_function = staticmethod(scf_input) 1608 input_required = False 1609 1610 1611class ScfForPhononsFactory(InputFactory): 1612 factory_function = staticmethod(scf_for_phonons) 1613 input_required = False 1614 1615 1616class PhononsFromGsFactory(InputFactory): 1617 factory_function = staticmethod(phonons_from_gsinput) 1618 1619 1620class PiezoElasticFactory(InputFactory): 1621 factory_function = staticmethod(scf_piezo_elastic_inputs) 1622 input_required = False 1623 1624 1625class PiezoElasticFromGsFactory(InputFactory): 1626 factory_function = staticmethod(piezo_elastic_inputs_from_gsinput) 1627