MICCoM School 2017 Ex#6 : Interface with ASE¶
You can try different systems, here are some examples.
We will use the ASE (Atomic Simulation Environment) Python package to build the input files.
6.1 : Benzene molecule¶
[1]:
# We download the pseudopotentials
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/C_ONCV_PBE-1.2.upf
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/H_ONCV_PBE-1.2.upf
[2]:
# We build the input for DFT
from ase import Atom, Atoms
from ase.cell import Cell
from ase.io import write
cell = Cell([[10,0,0],[0,10,0],[0,0,10]]) # Angstroms
positions = []
positions.append(Atom('C',[1.2115,0.6995,0.0000])) #Angstroms
positions.append(Atom('C',[1.2115,-0.6995,0.0000]))
positions.append(Atom('C',[0.0000,1.3990,0.0000]))
positions.append(Atom('C',[0.0000,-1.3990,0.0000]))
positions.append(Atom('C',[-1.2115,-0.6995,0.0000]))
positions.append(Atom('C',[-1.2115,0.6995,0.0000]))
positions.append(Atom('H',[0.0000,2.5000,0.0000]))
positions.append(Atom('H',[2.1651,1.2500,0.0000]))
positions.append(Atom('H',[2.1651,-1.2500,0.0000]))
positions.append(Atom('H',[0.0000,-2.5000,0.0000]))
positions.append(Atom('H',[-2.1651,-1.2500,0.0000]))
positions.append(Atom('H',[-2.1651,1.2500,0.0000]))
atoms = Atoms(positions,cell=cell)
pseudopotentials = {}
pseudopotentials['C'] = 'C_ONCV_PBE-1.2.upf'
pseudopotentials['H'] = 'H_ONCV_PBE-1.2.upf'
inp_data = {}
inp_data['calculation'] = 'scf'
inp_data['restart_mode'] = 'from_scratch'
inp_data['prefix'] = 'benzene'
inp_data['outdir'] = './'
inp_data['pseudo_dir'] = './'
inp_data['verbosity'] = 'high'
inp_data['wf_collect'] = True
inp_data['ecutwfc'] = 50
inp_data['nspin'] = 1
inp_data['nbnd'] = 40
inp_data['assume_isolated'] = 'mp'
inp_data['diago_full_acc'] = True
write('benzene.pwi',atoms,pseudopotentials=pseudopotentials,input_data=inp_data)
[3]:
!cat benzene.pwi
&CONTROL
calculation = 'scf'
verbosity = 'high'
restart_mode = 'from_scratch'
wf_collect = .true.
outdir = './'
prefix = 'benzene'
pseudo_dir = './'
/
&SYSTEM
nbnd = 40
ecutwfc = 50
nspin = 1
assume_isolated = 'mp'
ntyp = 2
nat = 12
ibrav = 0
/
&ELECTRONS
diago_full_acc = .true.
/
&IONS
/
&CELL
/
ATOMIC_SPECIES
C 12.011 C_ONCV_PBE-1.2.upf
H 1.008 H_ONCV_PBE-1.2.upf
K_POINTS gamma
CELL_PARAMETERS angstrom
10.00000000000000 0.00000000000000 0.00000000000000
0.00000000000000 10.00000000000000 0.00000000000000
0.00000000000000 0.00000000000000 10.00000000000000
ATOMIC_POSITIONS angstrom
C 1.2115000000 0.6995000000 0.0000000000
C 1.2115000000 -0.6995000000 0.0000000000
C 0.0000000000 1.3990000000 0.0000000000
C 0.0000000000 -1.3990000000 0.0000000000
C -1.2115000000 -0.6995000000 0.0000000000
C -1.2115000000 0.6995000000 0.0000000000
H 0.0000000000 2.5000000000 0.0000000000
H 2.1651000000 1.2500000000 0.0000000000
H 2.1651000000 -1.2500000000 0.0000000000
H 0.0000000000 -2.5000000000 0.0000000000
H -2.1651000000 -1.2500000000 0.0000000000
H -2.1651000000 1.2500000000 0.0000000000
[4]:
# We build the input for MBPT
import yaml
west = {}
west['input_west'] = {}
west['input_west']['qe_prefix'] = 'benzene'
west['input_west']['west_prefix'] = 'benzene'
west['input_west']['outdir'] = './'
west['wstat_control'] = {}
west['wstat_control']['wstat_calculation'] = 'S'
west['wstat_control']['n_pdep_eigen'] = 50
west['wfreq_control'] = {}
west['wfreq_control']['wfreq_calculation'] = 'XWGQ'
west['wfreq_control']['n_pdep_eigen_to_use'] = 50
west['wfreq_control']['qp_bandrange'] = [1,5]
west['wfreq_control']['n_refreq'] = 300
west['wfreq_control']['ecut_refreq'] = 2.0
with open('west.in', 'w') as file:
yaml.dump(west, file, sort_keys=False)
[5]:
!cat west.in
input_west:
qe_prefix: benzene
west_prefix: benzene
outdir: ./
wstat_control:
wstat_calculation: S
n_pdep_eigen: 50
wfreq_control:
wfreq_calculation: XWGQ
n_pdep_eigen_to_use: 50
qp_bandrange:
- 1
- 5
n_refreq: 300
ecut_refreq: 2.0
[ ]:
!mpirun -n 8 pw.x -i benzene.pwi > benzene.pwo
!mpirun -n 8 wstat.x -i wstat.in > wstat.out
!mpirun -n 8 wfreq.x -i wfreq.in > wfreq.out
6.2 : H2O¶
[6]:
# We download the pseudopotentials
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/O_ONCV_PBE-1.2.upf
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/H_ONCV_PBE-1.2.upf
!wget -N -q http://west-code.org/doc/training/water/H2O.cif
[7]:
# We build the input for DFT, starting from a CIF file
from ase.io import read, write
atoms = read('H2O.cif')
pseudopotentials = {}
pseudopotentials['O'] = 'O_ONCV_PBE-1.2.upf'
pseudopotentials['H'] = 'H_ONCV_PBE-1.2.upf'
inp_data = {}
inp_data['calculation'] = 'scf'
inp_data['restart_mode'] = 'from_scratch'
inp_data['prefix'] = 'h2o'
inp_data['outdir'] = './'
inp_data['pseudo_dir'] = './'
inp_data['verbosity'] = 'high'
inp_data['wf_collect'] = True
inp_data['ecutwfc'] = 60
inp_data['nspin'] = 1
inp_data['nbnd'] = 120
inp_data['diago_full_acc'] = True
write('h2o.pwi',atoms,pseudopotentials=pseudopotentials,input_data=inp_data)
[8]:
!cat h2o.pwi
&CONTROL
calculation = 'scf'
verbosity = 'high'
restart_mode = 'from_scratch'
wf_collect = .true.
outdir = './'
prefix = 'h2o'
pseudo_dir = './'
/
&SYSTEM
nbnd = 120
ecutwfc = 60
nspin = 1
ntyp = 2
nat = 36
ibrav = 0
/
&ELECTRONS
diago_full_acc = .true.
/
&IONS
/
&CELL
/
ATOMIC_SPECIES
H 1.008 H_ONCV_PBE-1.2.upf
O 15.999 O_ONCV_PBE-1.2.upf
K_POINTS gamma
CELL_PARAMETERS angstrom
7.60356630000000 0.00000000000000 0.00000000000000
-3.80178315000000 6.58488157515925 0.00000000000000
0.00000000000000 0.00000000000000 7.14296200000000
ATOMIC_POSITIONS angstrom
H 2.5154346141 0.0000000000 1.4030777397
H -1.2577173070 2.1784302773 1.4030777397
H 2.5440658430 4.4064512978 1.4030777397
H 5.0881316859 0.0000000000 4.9745587397
H 5.0595004570 4.4064512978 4.9745587397
H 1.2577173070 2.1784302773 4.9745587397
H 3.4394504051 0.0000000000 0.0760439735
H 5.8838410974 2.9786514259 0.0760439735
H 2.0820579474 3.6062301493 0.0760439735
H 4.1641158949 0.0000000000 3.6475249735
H -2.0820579474 3.6062301493 3.6475249735
H 1.7197252026 2.9786514259 3.6475249735
H 1.7103728160 5.7785299019 6.9918669248
H -2.0577569496 5.1768429328 6.9918669248
H 0.3473841335 2.2143903156 6.9918669248
H 5.5121559660 0.8063516733 6.9918669248
H 1.7440262004 1.4080386423 6.9918669248
H 4.1491672835 4.3704912595 6.9918669248
H 2.0914103340 0.8063516733 3.4203859248
H 5.8595400996 1.4080386423 3.4203859248
H 3.4543990165 4.3704912595 3.4203859248
H -1.7103728160 5.7785299019 3.4203859248
H 2.0577569496 5.1768429328 3.4203859248
H -0.3473841335 2.2143903156 3.4203859248
O 2.4836441034 0.0000000000 0.3975629790
O -1.2418220517 2.1508988875 0.3975629790
O 2.5599610983 4.4339826877 0.3975629790
O 5.1199221966 0.0000000000 3.9690439790
O 5.0436052017 4.4339826877 3.9690439790
O 1.2418220517 2.1508988875 3.9690439790
O 5.0229539156 0.0000000000 6.6399189011
O -2.5114769578 4.3500056930 6.6399189011
O 1.2903061922 2.2348758822 6.6399189011
O 2.5806123844 0.0000000000 3.0684379011
O -1.2903061922 2.2348758822 3.0684379011
O 2.5114769578 4.3500056930 3.0684379011
[9]:
# We build the input for MBPT
import yaml
west = {}
west['input_west'] = {}
west['input_west']['qe_prefix'] = 'h2o'
west['input_west']['west_prefix'] = 'h2o'
west['input_west']['outdir'] = './'
west['wstat_control'] = {}
west['wstat_control']['wstat_calculation'] = 'S'
west['wstat_control']['n_pdep_eigen'] = 200
west['wfreq_control'] = {}
west['wfreq_control']['wfreq_calculation'] = 'XWGQ'
west['wfreq_control']['macropol_calculation'] = 'C'
# Note that for a solid we must set macropol_calculation to 'C', the default 'N' is good for molecules
west['wfreq_control']['n_pdep_eigen_to_use'] = 200
west['wfreq_control']['qp_bandrange'] = [20,30]
west['wfreq_control']['n_refreq'] = 300
west['wfreq_control']['ecut_refreq'] = 2.0
with open('west.in', 'w') as file:
yaml.dump(west, file, sort_keys=False)
[10]:
!cat west.in
input_west:
qe_prefix: h2o
west_prefix: h2o
outdir: ./
wstat_control:
wstat_calculation: S
n_pdep_eigen: 200
wfreq_control:
wfreq_calculation: XWGQ
macropol_calculation: C
n_pdep_eigen_to_use: 200
qp_bandrange:
- 20
- 30
n_refreq: 300
ecut_refreq: 2.0
[ ]:
!mpirun -n 8 pw.x -i h2o.pwi > h2o.pwo
!mpirun -n 8 wstat.x -i wstat.in > wstat.out
!mpirun -n 8 wfreq.x -i wfreq.in > wfreq.out
6.3 : SiC solid¶
[11]:
# We download the pseudopotentials
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/C_ONCV_PBE-1.2.upf
!wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/Si_ONCV_PBE-1.2.upf
[12]:
# We build the input for DFT
from ase.build import bulk
atoms = bulk('SiC', 'zincblende', a=4.36, cubic=True) # Angstrom
pseudopotentials = {}
pseudopotentials['C'] = 'C_ONCV_PBE-1.2.upf'
pseudopotentials['Si'] = 'Si_ONCV_PBE-1.2.upf'
inp_data = {}
inp_data['calculation'] = 'scf'
inp_data['restart_mode'] = 'from_scratch'
inp_data['prefix'] = 'sic'
inp_data['outdir'] = './'
inp_data['pseudo_dir'] = './'
inp_data['verbosity'] = 'high'
inp_data['wf_collect'] = True
inp_data['ecutwfc'] = 60
inp_data['nspin'] = 1
inp_data['nbnd'] = 100
inp_data['diago_full_acc'] = True
write('sic_scf.pwi',atoms,pseudopotentials=pseudopotentials,input_data=inp_data,kpts=[8,8,8])
inp_data['calculation'] = 'nscf'
inp_data['nosym'] = True
inp_data['noinv'] = True
# Note that for systems with k-points we run a nscf calculation making sure that symmetries are turned off
write('sic_nscf.pwi',atoms,pseudopotentials=pseudopotentials,input_data=inp_data,kpts=[2,2,2])
[13]:
!cat sic_scf.pwi
&CONTROL
calculation = 'scf'
verbosity = 'high'
restart_mode = 'from_scratch'
wf_collect = .true.
outdir = './'
prefix = 'sic'
pseudo_dir = './'
/
&SYSTEM
nbnd = 100
ecutwfc = 60
nspin = 1
ntyp = 2
nat = 8
ibrav = 0
/
&ELECTRONS
diago_full_acc = .true.
/
&IONS
/
&CELL
/
ATOMIC_SPECIES
Si 28.085 Si_ONCV_PBE-1.2.upf
C 12.011 C_ONCV_PBE-1.2.upf
K_POINTS automatic
8 8 8 0 0 0
CELL_PARAMETERS angstrom
4.36000000000000 0.00000000000000 0.00000000000000
0.00000000000000 4.36000000000000 0.00000000000000
0.00000000000000 0.00000000000000 4.36000000000000
ATOMIC_POSITIONS angstrom
Si 0.0000000000 0.0000000000 0.0000000000
C 1.0900000000 1.0900000000 1.0900000000
Si 0.0000000000 2.1800000000 2.1800000000
C 1.0900000000 3.2700000000 3.2700000000
Si 2.1800000000 0.0000000000 2.1800000000
C 3.2700000000 1.0900000000 3.2700000000
Si 2.1800000000 2.1800000000 0.0000000000
C 3.2700000000 3.2700000000 1.0900000000
[14]:
!cat sic_nscf.pwi
&CONTROL
calculation = 'nscf'
verbosity = 'high'
restart_mode = 'from_scratch'
wf_collect = .true.
outdir = './'
prefix = 'sic'
pseudo_dir = './'
/
&SYSTEM
nbnd = 100
ecutwfc = 60
nosym = .true.
noinv = .true.
nspin = 1
ntyp = 2
nat = 8
ibrav = 0
/
&ELECTRONS
diago_full_acc = .true.
/
&IONS
/
&CELL
/
ATOMIC_SPECIES
Si 28.085 Si_ONCV_PBE-1.2.upf
C 12.011 C_ONCV_PBE-1.2.upf
K_POINTS automatic
2 2 2 0 0 0
CELL_PARAMETERS angstrom
4.36000000000000 0.00000000000000 0.00000000000000
0.00000000000000 4.36000000000000 0.00000000000000
0.00000000000000 0.00000000000000 4.36000000000000
ATOMIC_POSITIONS angstrom
Si 0.0000000000 0.0000000000 0.0000000000
C 1.0900000000 1.0900000000 1.0900000000
Si 0.0000000000 2.1800000000 2.1800000000
C 1.0900000000 3.2700000000 3.2700000000
Si 2.1800000000 0.0000000000 2.1800000000
C 3.2700000000 1.0900000000 3.2700000000
Si 2.1800000000 2.1800000000 0.0000000000
C 3.2700000000 3.2700000000 1.0900000000
[15]:
# We build the input for MBPT
import yaml
west = {}
west['input_west'] = {}
west['input_west']['qe_prefix'] = 'sic'
west['input_west']['west_prefix'] = 'sic'
west['input_west']['outdir'] = './'
west['wstat_control'] = {}
west['wstat_control']['wstat_calculation'] = 'S'
west['wstat_control']['n_pdep_eigen'] = 200
west['wfreq_control'] = {}
west['wfreq_control']['wfreq_calculation'] = 'XWGQ'
west['wfreq_control']['macropol_calculation'] = 'C'
west['wfreq_control']['n_pdep_eigen_to_use'] = 200
west['wfreq_control']['qp_bandrange'] = [20,30]
west['wfreq_control']['n_refreq'] = 300
west['wfreq_control']['ecut_refreq'] = 2.0
with open('west.in', 'w') as file:
yaml.dump(west, file, sort_keys=False)
[16]:
!cat west.in
input_west:
qe_prefix: sic
west_prefix: sic
outdir: ./
wstat_control:
wstat_calculation: S
n_pdep_eigen: 200
wfreq_control:
wfreq_calculation: XWGQ
macropol_calculation: C
n_pdep_eigen_to_use: 200
qp_bandrange:
- 20
- 30
n_refreq: 300
ecut_refreq: 2.0
[ ]:
!mpirun -n 8 pw.x -i sic_scf.pwi > sic_scf.pwo
!mpirun -n 8 pw.x -i sic_nscf.pwi > sic_nscf.pwo
!mpirun -n 8 wstat.x -i wstat.in > wstat.out
!mpirun -n 8 wfreq.x -i wfreq.in > wfreq.out