[1]:
# useful to autoreload the module without restarting the kernel
%load_ext autoreload
%autoreload 2
[2]:
from mppi import InputFiles as I
Tutorial for the PwInput class¶
This tutorial describes the main features and the usage of the PwInput class that enables to create and manage the input files for the QuantumESPRESSO computations
Setting up an input from scratch¶
Create an empty input object from scratch. If not parameters are provided some keys are init from the default data member of the class
[3]:
inp = I.PwInput()
inp
[3]:
{'control': {'calculation': "'scf'",
'verbosity': "'high'",
'prefix': "'pwscf'",
'outdir': "'./'"},
'system': {'force_symmorphic': '.false.'},
'electrons': {'diago_full_acc': '.false.'},
'ions': {},
'cell': {},
'atomic_species': {},
'atomic_positions': {},
'kpoints': {},
'cell_parameters': {}}
The input object can be initialized passing arguments in the constuctor both using the kwargs syntax and/or passing a dictionary. These operations override the default dictionary
[4]:
inp = I.PwInput(control = {'calculation' : 'nscf'})
inp
[4]:
{'control': {'calculation': 'nscf'},
'system': {'force_symmorphic': '.false.'},
'electrons': {'diago_full_acc': '.false.'},
'ions': {},
'cell': {},
'atomic_species': {},
'atomic_positions': {},
'kpoints': {},
'cell_parameters': {}}
[5]:
imp_dict = {'system' : {'ntyp' : 1}}
[6]:
inp = I.PwInput(control = {'calculation' : 'scf'}, **imp_dict)
inp
[6]:
{'control': {'calculation': 'scf'},
'system': {'ntyp': 1},
'electrons': {'diago_full_acc': '.false.'},
'ions': {},
'cell': {},
'atomic_species': {},
'atomic_positions': {},
'kpoints': {},
'cell_parameters': {}}
The attribute of the object can be updated using the standard procedure for python dictionaries, for instance
[7]:
inp['control'].update({'verbosity' : "'high'"})
inp
[7]:
{'control': {'calculation': 'scf', 'verbosity': "'high'"},
'system': {'ntyp': 1},
'electrons': {'diago_full_acc': '.false.'},
'ions': {},
'cell': {},
'atomic_species': {},
'atomic_positions': {},
'kpoints': {},
'cell_parameters': {}}
However it is more powerful to use the specific methods provided in the class (and other can be defined). For instance
[8]:
inp = I.PwInput()
inp.set_scf()
inp.set_prefix('si_scf_pref')
inp.set_pseudo_dir(pseudo_dir='pseudos')
# specific methods for the cell have still to be added
inp.set_lattice(ibrav=2,celldm1=10.3)
# Set the atoms type and positions
inp.add_atom('Si','Si.pbe-mt_fhi.UPF')
# add other atoms if needed, then set the number of atoms in the cell.
# This method sets the ntyp variable equal to the number of atoms added
# with the add_atom method
inp.set_atoms_number(2)
inp.set_atomic_positions([['Si',[0.,0.,0.,]],['Si',[0.25,0.25,0.25]]])
# Set the sampling of kpoints
inp.set_kpoints(type='automatic',points=[4,4,4])
inp
[8]:
{'control': {'calculation': "'scf'",
'verbosity': "'high'",
'prefix': "'si_scf_pref'",
'outdir': "'./'",
'pseudo_dir': "'pseudos'"},
'system': {'force_symmorphic': '.false.',
'ibrav': 2,
'celldm(1)': 10.3,
'ntyp': '1',
'nat': '2'},
'electrons': {'diago_full_acc': '.false.', 'conv_thr': 1e-08},
'ions': {},
'cell': {},
'atomic_species': {'Si': ['1.0', 'Si.pbe-mt_fhi.UPF']},
'atomic_positions': {'type': 'alat',
'values': [['Si', [0.0, 0.0, 0.0]], ['Si', [0.25, 0.25, 0.25]]]},
'kpoints': {'type': 'automatic', 'values': ([4, 4, 4], [0.0, 0.0, 0.0])},
'cell_parameters': {}}
The usage of the methods is described in the documentation and can be easily accessed as follows
[9]:
inp.set_pseudo_dir?
Signature: inp.set_pseudo_dir(pseudo_dir='pseudos', abs_path=False)
Docstring:
Set the position of the folder with the pseudo-potentials.
If `abs_path` is True the path is converted in a absolute path. In this way it
is possible to provide a relative path (expressed from the root of the folder where the notebook
is located) and the pseudo location can be found from an arbitrary folder.
Args:
pseudo_dir (:py:class:'string') : (relative) path of the folder with the pseduopotentials
Note:
If the folder tree contains blank spaces, QuantumESPRESSO cannot be able to find the pseudo, in this
cas it is safer to provide a relative path (expressed from the folder where the input file is
written)
File: ~/Applications/MPPI/mppi/InputFiles/PwInput.py
Type: method
[10]:
inp.set_kpoints?
Signature:
inp.set_kpoints(
type='automatic',
points=[1.0, 1.0, 1.0],
shift=[0.0, 0.0, 0.0],
klist=[],
)
Docstring:
Define the sampling of the Brillouin zone.
Args:
type (:py:class:`string`) : type of sampling (automatic, tpiba, tpiba_b,...)
points (:py:class:`list`) : number of kpoints in the x,y,z directions. Used only if
the type variable is set to `automatic`
shift (:py:class:`list`) : shifts in the x,y,z directions. Used only if the
type varible is set to `automatic`
klist(list) : list with the structure:
[[k1x,k1y,k1z,w1],[k2x,k2y,k2z,w2],....]
Used if type variable is not se to `automatic`
File: ~/Applications/MPPI/mppi/InputFiles/PwInput.py
Type: method
Once completed the input object can be converted to string with the correct format of the pw input files
[11]:
inp_tostring = inp.convert_string()
print(inp_tostring)
&control
calculation = 'scf'
outdir = './'
prefix = 'si_scf_pref'
pseudo_dir = 'pseudos'
verbosity = 'high'
/&end
&system
celldm(1) = 10.3
force_symmorphic = .false.
ibrav = 2
nat = 2
ntyp = 1
/&end
&electrons
conv_thr = 1e-08
diago_full_acc = .false.
/&end
ATOMIC_SPECIES
Si 1.0 Si.pbe-mt_fhi.UPF
ATOMIC_POSITIONS { alat }
Si 0.0000000000 0.0000000000 0.0000000000
Si 0.2500000000 0.2500000000 0.2500000000
K_POINTS { automatic }
4 4 4 0 0 0
and finally it can be written on file
[33]:
inp.write('IO_files/pw_from-scratch.in')
Build and input starting from an existing file¶
Consider a second example in which the input is initialized from an existing input file
[3]:
inp = I.PwInput(file='IO_files/gaas_scf.in')
inp
[3]:
{'control': {'calculation': "'scf'",
'verbosity': "'high'",
'prefix': "'ecut_40-k_4'",
'outdir': "'./'",
'pseudo_dir': "'../pseudos'"},
'system': {'force_symmorphic': '.false.',
'occupations': "'fixed'",
'ibrav': 2,
'celldm(1)': 10.677,
'ntyp': 2,
'nat': 2,
'ecutwfc': 40},
'electrons': {'diago_full_acc': '.false.', 'conv_thr': 1e-08},
'ions': {},
'cell': {},
'atomic_species': {'Ga': [69.72, 'Ga_hamlu.fhi.UPF'],
'As': [74.92, 'As_hamlu.fhi.UPF']},
'atomic_positions': {'type': 'alat',
'values': [['Ga', [0.0, 0.0, 0.0]], ['As', [0.25, 0.25, 0.25]]]},
'kpoints': {'type': 'automatic',
'values': ([4.0, 4.0, 4.0], [0.0, 0.0, 0.0])},
'cell_parameters': {},
'file': 'IO_files/gaas_scf.in'}
The variables can be modified, for instance we can modify the input to perform a nscf computation with 12 bands and a given k-path
[35]:
G = [0.,0.,0.]
X = [1.,0.,0.]
L = [0.5,0.5,0.5]
W = [1.0,0.5,0.]
[36]:
from mppi.Utilities import build_kpath
kpath = build_kpath(G,X,L,W,numstep=30)
kpath
[36]:
[[0.0, 0.0, 0.0, 30],
[1.0, 0.0, 0.0, 30],
[0.5, 0.5, 0.5, 30],
[1.0, 0.5, 0.0, 0]]
[38]:
inp.set_nscf(12)
inp.set_kpoints(type='tpiba_b',klist=kpath)
[39]:
inp_tostring = inp.convert_string()
print(inp_tostring)
&control
calculation = 'nscf'
outdir = './'
prefix = 'ecut_40-k_4'
pseudo_dir = '../pseudos'
verbosity = 'high'
/&end
&system
celldm(1) = 10.677
ecutwfc = 40
force_symmorphic = .false.
ibrav = 2
nat = 2
nbnd = 12
ntyp = 2
occupations = 'fixed'
/&end
&electrons
conv_thr = 1e-08
diago_full_acc = .false.
/&end
ATOMIC_SPECIES
Ga 69.72 Ga_hamlu.fhi.UPF
As 74.92 As_hamlu.fhi.UPF
ATOMIC_POSITIONS { alat }
Ga 0.0000000000 0.0000000000 0.0000000000
As 0.2500000000 0.2500000000 0.2500000000
K_POINTS { tpiba_b }
4
0.00000000 0.00000000 0.00000000 30.00000000
1.00000000 0.00000000 0.00000000 30.00000000
0.50000000 0.50000000 0.50000000 30.00000000
1.00000000 0.50000000 0.00000000 0.00000000
Finally a last example in which the lattice is specified through the ‘cell_parameters’ key
[40]:
inp = I.PwInput(file='IO_files/graphene_nscf.in')
inp
[40]:
{'control': {'calculation': "'nscf'",
'verbosity': "'high'",
'prefix': "'ecut",
'outdir': "'./'",
'pseudo_dir': "'../pseudos'"},
'system': {'force_symmorphic': '.false.',
'occupations': "'smearing'",
'smearing': "'fermi-dirac'",
'degauss': '0.0036749326',
'ibrav': '0',
'ntyp': '1',
'nat': '2',
'ecutwfc': '100',
'nbnd': '8'},
'electrons': {'diago_full_acc': '.false.', 'conv_thr': '1e-08'},
'ions': {},
'cell': {},
'atomic_species': {'C': ['12.011', 'C_pbe-20082014.UPF']},
'atomic_positions': {'type': 'angstrom',
'values': [['C', [0.0, 0.0, 0.0]], ['C', [0.0, 1.42, 0.0]]]},
'kpoints': {'type': 'tpiba_b',
'values': [[0.0, 0.0, 0.0, 40.0],
[0.5, 0.28867513, 0.0, 40.0],
[0.66666667, 0.0, 0.0, 40.0],
[0.0, 0.0, 0.0, 0.0]]},
'cell_parameters': {'type': 'angstrom',
'values': [[2.4595121467, 0.0, 0.0],
[1.2297560734, 2.13, 0.0],
[0.0, 0.0, 10.0]]},
'file': 'IO_files/graphene_nscf.in'}
[41]:
print(inp.convert_string())
&control
calculation = 'nscf'
outdir = './'
prefix = 'ecut
pseudo_dir = '../pseudos'
verbosity = 'high'
/&end
&system
degauss = 0.0036749326
ecutwfc = 100
force_symmorphic = .false.
ibrav = 0
nat = 2
nbnd = 8
ntyp = 1
occupations = 'smearing'
smearing = 'fermi-dirac'
/&end
&electrons
conv_thr = 1e-08
diago_full_acc = .false.
/&end
ATOMIC_SPECIES
C 12.011 C_pbe-20082014.UPF
ATOMIC_POSITIONS { angstrom }
C 0.0000000000 0.0000000000 0.0000000000
C 0.0000000000 1.4200000000 0.0000000000
K_POINTS { tpiba_b }
4
0.00000000 0.00000000 0.00000000 40.00000000
0.50000000 0.28867513 0.00000000 40.00000000
0.66666667 0.00000000 0.00000000 40.00000000
0.00000000 0.00000000 0.00000000 0.00000000
CELL_PARAMETERS angstrom
2.4595121467 0.0000000000 0.0000000000
1.2297560734 2.1300000000 0.0000000000
0.0000000000 0.0000000000 10.0000000000
[ ]: