xaxis and zaxis are not orthogonal !

[biegert]

Dear All,

I am trying to set up a calculation that will ultimately give me the imaginary part of the e-ph self energy for graphite - analog to PRL 112,257402 (2014) Fig. 2.

I am familiar with QE and I worked through some examples of EPW. I then tried to adapt the diamond example to graphite, made an identical directory structure with phonons/ epw/, set up scf and ph for graphite, ran pp.py (btw, small mod necessary to run on python3), copied the list of q-points from graphite.dyn0, made the positive k point list, and ran epw/scf and epw/nscf. All runs fine, but I get stuck in epw.in.

I suspect something is wrong with the projections or some parameters as the program crashes in the Wannier90 step after finding the k-point neighbors:

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Error in routine set_u_matrix (1):
xaxis and zaxis are not orthogonal !
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

stopping ...

I suspected that some of my projections are wrong, so I switched to 'random'. The code goes further but then crashes with:

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Error in routine createkmap (1):
q-vec not commensurate
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

I attach the inputs and outputs. I am really stuck here and tried quite a few different resolutions, convergence criteria, checked that the pseudo potential is norm conserving, etc etc.
I would really appreciate any hint what to do here to get this running.

With best regards,

Jens Biegert


------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
epw.in with "given projections"

--
&inputepw
prefix = 'graphite'
amass(1) = 12.01078
outdir = './'
dvscf_dir = '../phonons/save'
iverbosity = 0

elph = .true.
epbwrite = .true.
epbread = .false.

epwwrite = .true.
epwread = .false.

nbndsub = 8
nbndskip = 4

wannierize = .true.
num_iter = 300
iprint = 2
dis_win_max = 12
dis_win_min = -3
dis_froz_max= 7
dis_froz_min= -3
proj(1) = 'C:pz;sp2'

elecselfen = .false.
phonselfen = .false.
a2f = .false.

fsthick = 2.0 ! eV
eptemp = 300 ! K (same as PRB 76, 165108)
degaussw = 0.1 ! eV

nkf1 = 8
nkf2 = 8
nkf3 = 8
nqf1 = 6
nqf2 = 6
nqf3 = 6

nk1 = 8
nk1 = 8
nk2 = 8
nk3 = 8
nq1 = 6
nq2 = 6
nq3 = 6
/
28 cartesian
0.000000000000000E+00 0.000000000000000E+00 0.000000000000000E+00
etc...

------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
epw.out for "given projections"

Program EPW v.5.1.0 starts on 8Sep2019 at 17:46: 9

This program is part of the open-source Quantum ESPRESSO suite
for quantum simulation of materials; please cite
"P. Giannozzi et al., J. Phys.:Condens. Matter 21 395502 (2009);
"P. Giannozzi et al., J. Phys.:Condens. Matter 29 465901 (2017);
URL http://www.quantum-espresso.org",
in publications or presentations arising from this work. More details at
http://www.quantum-espresso.org/quote

Parallel version (MPI & OpenMP), running on 128 processor cores
Number of MPI processes: 32
Threads/MPI process: 4

MPI processes distributed on 1 nodes
K-points division: npool = 32

Reading data from directory:
./graphite.save/

IMPORTANT: XC functional enforced from input :
Exchange-correlation = PZ ( 1 1 0 0 0 0)
Any further DFT definition will be discarded
Please, verify this is what you really want

file C.UPF: wavefunction(s) 3d renormalized

G-vector sticks info
--------------------
sticks: dense smooth PW G-vecs: dense smooth PW
Sum 475 475 199 22807 22807 5847


Possibly too few bands at point 1 0.00000 0.00000 0.00000

Possibly too few bands at point 2 0.00000 0.00000 0.04585

Possibly too few bands at point 3 0.00000 0.00000 0.09170

Possibly too few bands at point 4 0.00000 0.00000 0.13754

Possibly too few bands at point 5 0.00000 0.00000 0.18339

Possibly too few bands at point 6 0.00000 0.00000 0.22924

Possibly too few bands at point 7 0.00000 0.00000 0.27509

Possibly too few bands at point 8 0.00000 0.00000 0.32094

Possibly too few bands at point 9 0.00000 0.14434 0.00000

Possibly too few bands at point 10 0.00000 0.14434 0.04585

Possibly too few bands at point 11 0.00000 0.14434 0.09170

Possibly too few bands at point 12 0.00000 0.14434 0.13754

Possibly too few bands at point 13 0.00000 0.14434 0.18339

Possibly too few bands at point 14 0.00000 0.14434 0.22924

Possibly too few bands at point 15 0.00000 0.14434 0.27509

Possibly too few bands at point 16 0.00000 0.14434 0.32094

--

bravais-lattice index = 4
lattice parameter (a_0) = 4.6412 a.u.
unit-cell volume = 236.0475 (a.u.)^3
number of atoms/cell = 4
number of atomic types = 1
kinetic-energy cut-off = 80.0000 Ry
charge density cut-off = 320.0000 Ry
Exchange-correlation = PZ ( 1 1 0 0 0 0)

celldm(1)= 4.64117 celldm(2)= 0.00000 celldm(3)= 2.72638
celldm(4)= 0.00000 celldm(5)= 0.00000 celldm(6)= 0.00000

crystal axes: (cart. coord. in units of a_0)
a(1) = ( 1.0000 0.0000 0.0000 )
a(2) = ( -0.5000 0.8660 0.0000 )
a(3) = ( 0.0000 0.0000 2.7264 )

reciprocal axes: (cart. coord. in units 2 pi/a_0)
b(1) = ( 1.0000 0.5774 0.0000 )
b(2) = ( 0.0000 1.1547 0.0000 )
b(3) = ( 0.0000 0.0000 0.3668 )


Atoms inside the unit cell:

Cartesian axes

site n. atom mass positions (a_0 units)
1 C 12.0108 tau( 1) = ( 0.00000 0.00000 0.68160 )
2 C 12.0108 tau( 2) = ( 0.00000 0.00000 2.04479 )
3 C 12.0108 tau( 3) = ( 0.00000 0.57735 0.68160 )
4 C 12.0108 tau( 4) = ( 0.50000 0.28868 2.04479 )

2 Sym.Ops. (with q -> -q+G )

G cutoff = 174.6002 ( 22807 G-vectors) FFT grid: ( 27, 27, 80)
number of k points= 512 gaussian broad. (Ry)= 0.0200 ngauss = 1
cart. coord. in units 2pi/a_0
k( 1) = ( 0.0000000 0.0000000 0.0000000), wk = 0.0039062
k( 2) = ( 0.0000000 0.0000000 0.0458483), wk = 0.0039062
etc...

PseudoPot. # 1 for C read from file:
/opt/qe-6.4.1/pseudo/C.UPF
MD5 check sum: 4781e8ce5ee01a432381f51d0ea25c53
Pseudo is Norm-conserving, Zval = 4.0
Generated by new atomic code, or converted to UPF format
Using radial grid of 461 points, 2 beta functions with:
l(1) = 0
l(2) = 1
EPW : 2.79s CPU 2.43s WALL

EPW : 3.62s CPU 2.68s WALL

No wavefunction gauge setting applied
-------------------------------------------------------------------
Wannierization on 8 x 8 x 8 electronic grid
-------------------------------------------------------------------

Spin CASE ( default = unpolarized )

Initializing Wannier90


Initial Wannier projections

( 0.00000 0.00000 0.25000) : l = 0 mr = 1
( 0.00000 0.00000 0.25000) : l = 0 mr = 2
( 0.00000 0.00000 0.25000) : l = 0 mr = 3
( 0.00000 0.00000 0.25000) : l = 0 mr = 1
( 0.00000 0.00000 0.75000) : l = 0 mr = 1
( 0.00000 0.00000 0.75000) : l = 0 mr = 2
( 0.00000 0.00000 0.75000) : l = 0 mr = 3
( 0.00000 0.00000 0.75000) : l = 0 mr = 1

- Number of bands is ( 8)
- Number of total bands is ( 8)
- Number of excluded bands is ( 0)
- Number of wannier functions is ( 8)
- All guiding functions are given

Reading data about k-point neighbours

- All neighbours are found


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Error in routine set_u_matrix (1):
xaxis and zaxis are not orthogonal !
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

stopping ...

------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
epw.in with "random projections"

same as epw.in above but with
proj(1) = 'random'

------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
epw.out with "random projections"; the calculated Wannier centers do not make sense to me here.

same as above and then

Initial Wannier projections

( 0.78835 0.71458 0.81618) : l = 0 mr = 1
( 0.41054 0.98136 0.63028) : l = 0 mr = 1
( 0.40172 0.90805 0.39456) : l = 0 mr = 1
( 0.23740 0.79789 0.76079) : l = 0 mr = 1
( 0.56474 0.23536 0.30491) : l = 0 mr = 1
( 0.25089 0.62963 0.17798) : l = 0 mr = 1
( 0.55769 0.47664 0.98269) : l = 0 mr = 1
( 0.50964 0.18660 0.17143) : l = 0 mr = 1

- Number of bands is ( 8)
- Number of total bands is ( 8)
- Number of excluded bands is ( 0)
- Number of wannier functions is ( 8)
- All guiding functions are given

Reading data about k-point neighbours

- All neighbours are found

AMN
k points = 512 in 32 pools
1 of 16 on ionode
2 of 16 on ionode
3 of 16 on ionode
4 of 16 on ionode
5 of 16 on ionode
6 of 16 on ionode
7 of 16 on ionode
8 of 16 on ionode
9 of 16 on ionode
10 of 16 on ionode
11 of 16 on ionode
12 of 16 on ionode
13 of 16 on ionode
14 of 16 on ionode
15 of 16 on ionode
16 of 16 on ionode

AMN calculated

MMN
k points = 512 in 32 pools
1 of 16 on ionode
2 of 16 on ionode
3 of 16 on ionode
4 of 16 on ionode
5 of 16 on ionode
6 of 16 on ionode
7 of 16 on ionode
8 of 16 on ionode
9 of 16 on ionode
10 of 16 on ionode
11 of 16 on ionode
12 of 16 on ionode
13 of 16 on ionode
14 of 16 on ionode
15 of 16 on ionode
16 of 16 on ionode
MMN calculated

Running Wannier90

Wannier Function centers (cartesian, alat) and spreads (ang):

( 0.42616 0.56789 2.05054) : 1.92075
( -0.06504 1.05902 1.92261) : 1.94633
( -0.16375 0.68767 0.78751) : 1.37537
( -0.32848 0.44933 2.05039) : 1.96696
( 0.22356 0.67433 0.70294) : 0.78409
( -0.04258 0.35766 0.59651) : 1.03597
( 0.42854 0.18108 2.15558) : 1.89177
( 0.47768 -0.27536 0.63942) : 1.93636

-------------------------------------------------------------------
WANNIER : 468.97s CPU 142.98s WALL ( 1 calls)
-------------------------------------------------------------------

Dipole matrix elements calculated

Calculating kgmap

Progress kgmap: ########################################
kmaps : 111.66s CPU 115.88s WALL ( 1 calls)
Symmetries of Bravais lattice: 24
Symmetries of crystal: 24


===================================================================
irreducible q point # 1
===================================================================

Symmetries of small group of q: 24
in addition sym. q -> -q+G:

Number of q in the star = 1
List of q in the star:
1 0.000000000 0.000000000 0.000000000
Imposing acoustic sum rule on the dynamical matrix

q( 1 ) = ( 0.0000000 0.0000000 0.0000000 )


===================================================================
irreducible q point # 2
===================================================================

Symmetries of small group of q: 12

Number of q in the star = 2
List of q in the star:
1 0.000000000 0.000000000 0.061131023
2 0.000000000 0.000000000 -0.061131023

q( 2 ) = ( 0.0000000 0.0000000 0.0611310 )

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Error in routine createkmap (1):
q-vec not commensurate
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

stopping ...










------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
scf.in ; BTW I think that using nbnd here is wrong as it is only used in nscf, but I followed the given examples

&control
calculation = 'scf'
prefix = 'graphite'
restart_mode = 'from_scratch'
wf_collect = .false.
pseudo_dir = '/opt/qe-6.4.1/pseudo/'
outdir = './'
tprnfor = .true.
tstress = .true.
/
&system
ibrav = 4
a = 2.456
c = 6.696
nat = 4
ntyp = 1
ecutwfc = 80
occupations = 'smearing'
smearing = 'mp'
degauss = 0.02
nbnd = 8
/
&electrons
diagonalization = 'david'
mixing_beta = 0.7
conv_thr = 1.0d-10
/
ATOMIC_SPECIES
C 12.01078 C.UPF
ATOMIC_POSITIONS crystal
C 0.000000000000 0.000000000000 0.250000000000
C 0.000000000000 0.000000000000 0.750000000000
C 0.333333333333 0.666666666666 0.250000000000
C 0.666666666666 0.333333333333 0.750000000000
K_POINTS automatic
8 8 8 0 0 0

------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
nscf.in

&control
calculation = 'nscf'
prefix = 'graphite'
pseudo_dir = '/opt/qe-6.4.1/pseudo/'
outdir = './'
/
&system
ibrav = 4
a = 2.456
c = 6.696
nat = 4
ntyp = 1
ecutwfc = 80
occupations = 'smearing'
smearing = 'mp'
degauss = 0.02
nbnd = 8
nosym = .true.
/
&electrons
diagonalization = 'david'
mixing_beta = 0.7
conv_thr = 1.0d-10
/
ATOMIC_SPECIES
C 12.01078 C.UPF
ATOMIC_POSITIONS crystal
C 0.000000000000 0.000000000000 0.250000000000
C 0.000000000000 0.000000000000 0.750000000000
C 0.333333333333 0.666666666666 0.250000000000
C 0.666666666666 0.333333333333 0.750000000000
K_POINTS crystal
512
0.00000000 0.00000000 0.00000000 1.953125e-03
0.00000000 0.00000000 0.12500000 1.953125e-03

[roxana]

Hi,

Can you try first to get the wannierization using the Wannier90 code in order to single out the problem?

Best,
Roxana