Commit 4132a647 authored by Pietro Delugas's avatar Pietro Delugas

Merge branch 'master_641' into 'master'

Master 641

See merge request !396
parents 559c6f52 a04f648e
......@@ -148,8 +148,7 @@ MODULE cp_restart_new
CHARACTER(LEN=20) :: dft_name
CHARACTER(LEN=256) :: dirname
CHARACTER(LEN=320) :: filename, sourcefile
CHARACTER(LEN=4) :: cspin
INTEGER :: kunit, ik_eff
INTEGER :: ik_eff
INTEGER :: k1, k2, k3
INTEGER :: nk1, nk2, nk3
INTEGER :: j, i, iss, ig, nspin_wfc, iss_wfc
......@@ -646,11 +645,8 @@ MODULE cp_restart_new
COMPLEX(DP), INTENT(INOUT) :: cm2(:,:) !
REAL(DP), INTENT(INOUT) :: wfc(:,:) ! BS
!
CHARACTER(LEN=256) :: dirname, kdirname, filename
CHARACTER(LEN=5) :: kindex
CHARACTER(LEN=4) :: cspin
CHARACTER(LEN=256) :: dirname, filename
INTEGER :: strlen
INTEGER :: kunit
INTEGER :: k1, k2, k3
INTEGER :: nk1, nk2, nk3
INTEGER :: i, j, iss, ig, nspin_wfc, ierr, ik
......
......@@ -702,13 +702,17 @@ makov_payne.o : ../../UtilXlib/mp.o
makov_payne.o : ../../UtilXlib/parallel_include.o
makov_payne.o : ions_positions.o
makov_payne.o : mainvar.o
manycp.o : ../../LAXlib/mp_diag.o
manycp.o : ../../Modules/check_stop.o
manycp.o : ../../Modules/command_line_options.o
manycp.o : ../../Modules/environment.o
manycp.o : ../../Modules/input_parameters.o
manycp.o : ../../Modules/io_global.o
manycp.o : ../../Modules/mp_bands.o
manycp.o : ../../Modules/mp_global.o
manycp.o : ../../Modules/mp_images.o
manycp.o : ../../Modules/mp_pools.o
manycp.o : ../../Modules/mp_world.o
manycp.o : ../../Modules/read_input.o
manycp.o : input.o
metaxc.o : ../../Modules/funct.o
......
No preview for this file type
......@@ -50,6 +50,7 @@ calculation to define paths in the BZ. This feature is available with
the option \texttt{tpiba\_b} or \texttt{crystal\_b} in a \texttt{'bands'}
calculation or with the option \texttt{q\_in\_band\_form} in the input of the
\texttt{matdyn.x} code.
BEWARE: you need to explicitly specify \texttt{ibrav} to use this feature.
Lines in reciprocal space are defined by giving the coordinates of the
starting and ending points and the number of points of each line.
The coordinates of the starting and ending points can be
......
No preview for this file type
\documentclass[12pt,a4paper]{article}
\def\version{6.4}
\def\version{6.4.1}
\def\QE{{\sc Quantum ESPRESSO}}
\def\qe{QE}
\textwidth = 17cm
......@@ -1194,8 +1194,9 @@ may save hours of searching into the code for a piece of missing information.
complain if the latter \& is missing, others do.
% Another example: empty strings are nonstandard,
% use \texttt{empty='~'}, not \texttt{empty=''}.
\item do not (yet) use F2008 syntax. Stick to F2003 at most (for now).
\qe\ must work even if you do not have the latest and the greatest compiler.
\item try to stick to F2003 standard: \qe\ must work even if you do not have
the latest and the greatest compiler. Use F2008 syntax only if really
useful, and after verifying that it doesn't break too many compilers.
\item use "dp" (defined in module ''kinds'') to define the type of real and
complex variables
\item all constants should be defined to be of kind "dp". Preferred syntax:
......@@ -1210,7 +1211,7 @@ double precision complex number).
\item Do not use automatic arrays (e.g. \texttt{REAL(dp) :: A(N)} with
\texttt{N} defined at run time) unless you are sure that the array is
small in all cases: large arrays may easily exceed the stack size,
or the memory size,
or the memory size.
\item Do not use pointers unless you have a good reason to:
pointers may hinder optimization. Allocatable arrays should be used instead.
\item If you use pointers, nullify them before performing tests on their
......@@ -1224,9 +1225,12 @@ array sections. Passing an array section to a routine may look elegant
but it may turn out to be inefficient: a copy will be silently done
if the section is not contiguous in memory (or if the compiler
decides it is the right thing to do), increasing the memory footprint.
\item Do not pass unallocated arrays as arguments, even in those cases where
they are not actually used inside the subroutine: some compilers don't
like it.
\item Do not pass unallocated arrays or pointers as non-optional arguments,
even in those cases where they are not actually used inside the subroutine:
some compilers don't like it. Also note that if passed as optional argument
--provided the argument has not the pointer or allocatable attribute--
unallocated arrays or pointers are interpreted as non present (this is a
F2008 feature, already used since v.6.4).
\item Do not use any construct that is susceptible to be flagged as
out-of-bounds error, even if no actual out-of-bound error takes place.
\item Always use IMPLICIT NONE and declare all local variables.
......
New in 6.4.1 branch :
* A warning is issued if the lattice parameter seems to be a conversion
factor instead of a true lattice parameter. Conversion should be achieved
with the appropriate options, not with dirty tricks. In the future this
will no longer be allowed
* A warning is issued if ibrav=0 is used for systems having symmetry. If not
properly done this may lead to strange problems with symmetry detection
and symmetrization. Lattice information should be used if available.
Problems fixed in 6.4.1 branch :
* Two bugs fixed in HP: 1) the code was not working correctly when fractional
translations were present, 2) there was a bug in the case when either there
is only one k point, or when k pools are used and some of the pools have
only one k point.
* Restart of ph.x with 2D boundary conditions has been fixed (see gitlab
issue #102)
* XML file correctly written if tetrahedra are used (see gitlab issue #103)
New in version 6.4:
* Experimental version of SCDM localization with k-points, activated like for
k=0 by specifying in &system namelist a value > 0 for "localization_thr".
......@@ -17,6 +35,8 @@ New in version 6.4:
* XDM now works also for USPP and norm-conserving PP
Problems fixed in version 6.4 (+ = in qe-6.3-backports as well) :
+ Codes reading scf data recomputed celldm parameters also if ibrav=0
This produced confusing output and had the potential to break some codes
+ index not correctly initialized in LSDA phonon with core corrections
+ GTH pseudopotentials in analytical form wrongly computed in some cases
+ projwfc.x not working with new xml format in noncolinear/spinorbit case
......@@ -55,6 +75,8 @@ Incompatible changes in version 6.4 version:
Known problems in version 6.4:
* Frequent "dexx is negative" errors with hybrid functionals
* restart of ph.x when using 2D boundary conditions fails see issue#102 on gitLab
* the band_structure element is printed incompletely if tetrahedra are used for sums in the IBZ, see issue #103.
New in 6.3 version:
* New implementation, using a more robust algorithm for the Wigner-Seitz
......
No preview for this file type
\documentclass[12pt,a4paper]{article}
\def\version{6.4}
\def\version{6.4.1}
\def\qe{{\sc Quantum ESPRESSO}}
\usepackage{html}
......@@ -407,7 +407,9 @@ may break them. Use \texttt{export LC\_ALL=C} (sh/bash) or
when running scripts (including installation scripts).
Second, you need C and Fortran compilers, compliant with C89 and
F2003 standards. For parallel
F2003 standards\footnote{since v.6.4 a standard 2008 feature is
used: if unallocated pointers are passed as optional arguments,
they are interpreted as not present}. For parallel
execution, you will also need MPI libraries and a parallel
(i.e. MPI-aware) compiler. For massively parallel machines, or
for simple multicore parallelization, an OpenMP-aware compiler
......@@ -945,6 +947,7 @@ the last written output file to understand why.
\begin{itemize}
\item
Working Fortran and C compilers, compliant with F2003 and C89 standards
(see Sec.\ref{Sec:Installation})
respectively, are needed in order to compile \qe. Most recent Fortran
compilers will do the job.
......@@ -1109,6 +1112,10 @@ add preprocessing option \texttt{-Dzdotc=zdotc\_wrapper} to \texttt{DFLAGS}.
\paragraph{Linux PCs with Intel compiler (ifort)}
IMPORTANT NOTE: ifort versions earlier than v.15 miscompile the new
XML code in QE v.6.4 and later. Please install this patch:\\
\texttt{https://gitlab.com/QEF/q-e/wikis/Support/Patch-for-old-Intel-compilers}.
The Intel compiler ifort \texttt{http://software.intel.com/}
produces fast executables, at least on Intel CPUs, but not all versions
work as expected. In case of trouble, update your version
......@@ -1132,8 +1139,6 @@ The warning: {\em feupdateenv is not implemented and will always fail},
can be safely ignored. Warnings on ``bad preprocessing option'' when compiling
iotk and complains about ``recommended formats'' may also be ignored.
Versions v.12 and earlier of ifort are no longer supported by QE v.\version.
\paragraph{Linux PCs with MKL libraries}
On Intel CPUs it is very convenient to use Intel MKL libraries
(freely available at
......
This diff is collapsed.
......@@ -4,191 +4,11 @@
# 14/07/2015 - Creation of the script - Samuel Ponce
# 14/03/2018 - Automatically reads the number of q-points - Michael Waters
# 14/03/2018 - Detect if SOC is included in the calculation - Samuel Ponce
# 13/11/2018 - Write dyn files in xml format for SOC case - Shunhong Zhang (USTC)
#
import numpy as np
import os
from xml.dom import minidom
# Convert the dyn files to the xml form, for SOC case - Shunhong Zhang (USTC)
def dyn2xml(prefix):
ndyn=int(os.popen('head -2 {0}.dyn0|tail -1'.format(prefix)).read())
for idyn in range(1,ndyn+1):
print '{0}.dyn{1} to {0}.dyn_q{1}.xml'.format(prefix,idyn)
dynmat=dyn(prefix,idyn)
dynmat._write_xml()
def get_geom_info():
if os.path.isfile('ph.out')==False:
print 'cannot extract geometry info from ph.out'
return 1
else:
volm=float(os.popen('grep -a volume ph.out 2>/dev/null|tail -1').readline().split()[-2])
get_at=os.popen('grep -a -A 3 "crystal axes" ph.out 2>/dev/null|tail -3').readlines()
at=np.array([[float(item) for item in line.split()[3:6]] for line in get_at])
get_bg=os.popen('grep -a -A 3 "reciprocal axes" ph.out 2>/dev/null|tail -3').readlines()
bg=np.array([[float(item) for item in line.split()[3:6]] for line in get_bg])
return volm,at,bg
class dyn(object):
def __init__(self,prefix,idyn):
self._prefix=prefix
self._idyn=idyn
fil='{0}.dyn{1}'.format(prefix,idyn)
f=open(fil)
self._comment=f.readline()
f.readline()
line=f.readline().split()
self._ntype=int(line[0])
self._natom=int(line[1])
self._ibrav=int(line[2])
self._nspin=1
self._cell_dim=np.array([float(ii) for ii in line[3:]])
self._volm=0
self._at=np.zeros((3,3),float)
self._bg=np.zeros((3,3),float)
try: self._volm,self._at,self._bg = get_geom_info()
except: print 'warning: lattice info not found'
self._species=[];
self._mass=[]
for i in range(self._ntype):
line=f.readline().split()
self._species.append(line[1].strip("'"))
self._mass.append(float(line[-1])/911.4442) # normalize to atomic mass
self._atom_type=np.zeros(self._natom,int)
self._pos=np.zeros((self._natom,3),float)
for i in range(self._natom):
line=f.readline().split()
self._atom_type[i]=int(line[1])
for j in range(3): self._pos[i,j]=float(line[j+2])
self._nqpt=int(os.popen('grep -c "Dynamical Matrix" {0}'.format(fil)).read().split()[0])
self._qpt=[]
self._dynmat=np.zeros((self._nqpt,self._natom,self._natom,3,3,2),float)
f.readline()
for iqpt in range(self._nqpt):
f.readline();
f.readline()
line=f.readline().split()
self._qpt.append(np.array([float(item) for item in line[3:6]]))
f.readline()
for i in range(self._natom):
for j in range(self._natom):
f.readline()
data=np.fromfile(f,sep=' ',count=18,dtype=float).reshape(3,3,2)
self._dynmat[iqpt,i,j]=data
self._qpt=np.array(self._qpt)
for i in range(5): f.readline()
self._freq=np.zeros((self._natom*3,2),float)
self._disp=np.zeros((self._natom*3,self._natom,3,2),float)
for i in range(self._natom*3):
line=f.readline().split()
self._freq[i,0]=float(line[4])
self._freq[i,1]=float(line[7])
for j in range(self._natom):
line=f.readline().split()[1:-1]
data=np.array([float(item) for item in line]).reshape(3,2)
self._disp[i,j]=data
def _write_xml(self):
doc=minidom.Document()
root = doc.createElement('Root')
doc.appendChild(root)
geom_info=doc.createElement('GEOMETRY_INFO')
tags=('NUMBER_OF_TYPES','NUMBER_OF_ATOMS','BRAVAIS_LATTICE_INDEX','SPIN_COMPONENTS')
numbers=(self._ntype,self._natom,self._ibrav,self._nspin)
for i,(tag,num) in enumerate(zip(tags,numbers)):
inode=doc.createElement(tag)
inode.setAttribute('type','integer')
inode.setAttribute('size','1')
inode.text=num
inode.appendChild(doc.createTextNode(str(num)))
geom_info.appendChild(inode)
cell_dim=doc.createElement('CELL_DIMENSIONS')
cell_dim.setAttribute('type','real')
cell_dim.setAttribute('size','6')
for i in range(6):
cell_dim.appendChild(doc.createTextNode('{0:16.10f}'.format(self._cell_dim[i])))
geom_info.appendChild(cell_dim)
tags=['AT','BG']
for tag,lat in zip(tags,(self._at,self._bg)):
inode=doc.createElement(tag)
inode.setAttribute('type','real')
inode.setAttribute('size','9')
inode.setAttribute('columns','3')
for i in range(3):
text=' '.join(['{0:16.10f}'.format(item) for item in lat[i]])
inode.appendChild(doc.createTextNode(text))
geom_info.appendChild(inode)
volm=doc.createElement('UNIT_CELL_VOLUME_AU')
volm.setAttribute('type','real')
volm.setAttribute('size','1')
volm.appendChild(doc.createTextNode('{0:16.10f}'.format(self._volm)))
geom_info.appendChild(volm)
for itype in range(self._ntype):
nt=doc.createElement('TYPE_NAME.{0}'.format(itype+1))
nt.setAttribute('type','character')
nt.setAttribute('size','1')
nt.setAttribute('len','3')
nt.appendChild(doc.createTextNode('{0}'.format(self._species[itype])))
na=doc.createElement('MASS.{0}'.format(itype+1))
na.setAttribute('type','real')
na.setAttribute('size','1')
na.appendChild(doc.createTextNode('{0:16.10f}'.format(self._mass[itype])))
geom_info.appendChild(nt)
geom_info.appendChild(na)
for iat in range(self._natom):
at=doc.createElement('ATOM.{0}'.format(iat+1))
at.setAttribute('SPECIES','{0}'.format(self._species[self._atom_type[iat]-1]))
at.setAttribute('INDEX',str(iat+1))
pos=' '.join(['{0:16.10f}'.format(item) for item in self._pos[iat]])
at.setAttribute('TAU',pos)
geom_info.appendChild(at)
nqpt=doc.createElement('NUMBER_OF_Q')
nqpt.setAttribute('type','integer')
nqpt.setAttribute('size','1')
nqpt.appendChild(doc.createTextNode(str(self._nqpt)))
geom_info.appendChild(nqpt)
root.appendChild(geom_info)
for iqpt in range(self._nqpt):
dynmat=doc.createElement('DYNAMICAL_MAT_.{0}'.format(iqpt+1))
qpt=doc.createElement('Q_POINT')
qpt.setAttribute('type','real')
qpt.setAttribute('size','3')
qpt.setAttribute('columns','3')
tnode=doc.createTextNode(' '.join(['{0:16.10f}'.format(item) for item in self._qpt[iqpt]]))
qpt.appendChild(tnode)
dynmat.appendChild(qpt)
for iat in range(self._natom):
for jat in range(self._natom):
ph=doc.createElement('PHI.{0}.{1}'.format(iat+1,jat+1))
ph.setAttribute('type','complex')
ph.setAttribute('size','9')
ph.setAttribute('columns','3')
for i in range(3):
for j in range(3):
text='{0:16.10f} {1:16.10f}'.format(self._dynmat[iqpt,iat,jat,i,j,0],self._dynmat[iqpt,iat,jat,i,j,1])
ph.appendChild(doc.createTextNode(text))
dynmat.appendChild(ph)
root.appendChild(dynmat)
mode=doc.createElement('FREQUENCIES_THZ_CMM1')
for iomega in range(self._natom*3):
inode=doc.createElement('OMEGA.{0}'.format(iomega+1))
inode.setAttribute('type','real')
inode.setAttribute('size','2')
inode.setAttribute('columns','2')
inode.appendChild(doc.createTextNode('{0:16.10f} {1:16.10f}'.format(self._freq[iomega,0],self._freq[iomega,1])))
idisp=doc.createElement('DISPLACEMENT.{0}'.format(iomega+1))
idisp.setAttribute('tpye','complex')
idisp.setAttribute('size','3')
for iat in range(self._natom):
for j in range(3):
tnode=doc.createTextNode('{0:16.10f} {1:16.10f}'.format(self._disp[iomega,iat,j,0],self._disp[iomega,iat,j,1]))
idisp.appendChild(tnode)
mode.appendChild(inode)
mode.appendChild(idisp)
root.appendChild(mode)
fp = open('{0}.dyn_q{1}.xml'.format(self._prefix,self._idyn), 'w')
doc.writexml(fp, addindent=' ', newl='\n')
# Return the number of q-points in the IBZ
def get_nqpt(prefix):
fname = '_ph0/' +prefix+'.phsave/control_ph.xml'
......@@ -233,21 +53,6 @@ prefix = str(user_input)
# Test if SOC
SOC = hasSOC(prefix)
# If SOC detected, but dyn is not in XML and we want to convert it
if SOC=='true':
user_input = raw_input('Calculation with SOC detected. Do you want to convert dyn in XML format [y/n]?\n')
if str(user_input) == 'y':
dyn2xml(prefix)
os.system('mv {0}.dyn*.xml save'.format(prefix))
# If no SOC detected, do you want to convert into XML format
if SOC=='false':
user_input = raw_input('Calculation without SOC detected. Do you want to convert to xml anyway [y/n]?\n')
if str(user_input) == 'y':
SOC = 'true'
dyn2xml(prefix)
os.system('mv {0}.dyn*.xml save'.format(prefix))
# Test if seq. or parallel run
SEQ = isSEQ(prefix)
......
......@@ -10,7 +10,6 @@ authors: Samuel Poncé
author_description: The EPW project is mainly developed at the university of Oxford.
github: https://github.com/sponce24
email: samuel.pon@gmail.com
project_sourceforge: http://qeforge.qe-forge.org/gf/project/q-e/
predocmark: >
media_dir: ./media
page_dir: ./Ford
......
......@@ -645,8 +645,6 @@
! every processor has just a chunk of the array, I may need some
! communication)
!
! No ultrasoft now
!
! I use the rule : if not found then gmap = 0
! Note that the map will be used only up to npwx (small sphere),
! while the G-vectors lost in the process are on the surface of
......
......@@ -65,8 +65,6 @@
zi_all(:,:), &!
esigmar_all(:,:,:), &!
esigmai_all(:,:,:), &!
gammar_all(:,:,:), &! Real part of the Phonon self-energy (freq. dependent for spectral function)
gammai_all(:,:,:), &! Imaginary part of the Phonon self-energy (freq. dependent for spectral function)
jdos(:), &!
spectra(:,:,:,:,:,:), &! dipole absorption spectra, polarizations, nomega, nsmear, dme/vme, absorption/emission
sumr(:,:,:,:), &! to apply the ASR correction to every xq
......
......@@ -18,9 +18,6 @@
!! Compact formalism Dec 2006
!! Phonon irreducible zone Mar 2007
!!
!! No ultrasoft now
!! No spin polarization
!!
!! RM - add noncolin case
!-----------------------------------------------------------------------
!
......@@ -1475,6 +1472,13 @@
!
ENDDO
CLOSE(lambda_phself)
!
! SP - 03/2019
! \Gamma = 1/\tau = phonon lifetime
! \Gamma = - 2 * Im \Pi^R where \Pi^R is the retarted phonon self-energy.
! Im \Pi^R = pi*k-point weight*[f(E_k+q) - f(E_k)]*delta[E_k+q - E_k - w_q]
! Since gamma_all = pi*k-point weight*[f(E_k) - f(E_k+q)]*delta[E_k+q - E_k - w_q] we have
! \Gamma = 2 * gamma_all
OPEN(unit=linewidth_phself,file='linewidth.phself')
WRITE(linewidth_phself, '(a)') '# Phonon frequency and phonon lifetime in meV '
WRITE(linewidth_phself,'(a)') '# Q-point Mode Phonon freq (meV) Phonon linewidth (meV)'
......@@ -1482,7 +1486,7 @@
!
DO imode=1, nmodes
WRITE(linewidth_phself,'(i9,i6,E20.8,E22.10)') iqq,imode,&
ryd2mev*wf(imode,iqq),ryd2mev*REAL(gamma_all(imode,iqq,1))
ryd2mev*wf(imode,iqq), 2.0d0 * ryd2mev * REAL(gamma_all(imode,iqq,1))
ENDDO
!
ENDDO
......
......@@ -20,9 +20,6 @@
!! Compact formalism Dec 2006
!! Phonon irreducible zone Mar 2007
!!
!! No ultrasoft now
!! No spin polarization
!!
!! RM - add noncolin case
!-----------------------------------------------------------------------
!
......
......@@ -20,8 +20,6 @@
!! every processor has just a chunk of the array, I may need some
!! communication)
!!
!! No ultrasoft now
!!
!----------------------------------------------------------------------
USE kinds, ONLY : DP
USE constants_epw, ONLY : twopi, ci, cone
......
......@@ -2074,8 +2074,8 @@
USE io_epw, ONLY : iufilgap
USE io_files, ONLY : prefix
USE epwcom, ONLY : fsthick
USE eliashbergcom, ONLY : estemp, Agap, nkfs, nbndfs, ef0, ekfs
USE constants_epw, ONLY : kelvin2eV, zero
USE eliashbergcom, ONLY : estemp, Agap, nkfs, nbndfs, ef0, ekfs, w0g
USE constants_epw, ONLY : kelvin2eV, zero, eps5
!
IMPLICIT NONE
!
......@@ -2099,8 +2099,6 @@
!! Step size in nbin
REAL(DP) :: delta_max
!! Max value of superconducting gap
REAL(DP) :: sigma
!! Variable for smearing
REAL(DP) :: weight
!! Variable for weight
REAL(DP), ALLOCATABLE :: delta_k_bin(:)
......@@ -2110,8 +2108,8 @@
!
temp = estemp(itemp) / kelvin2eV
!
delta_max = 1.25d0 * maxval(Agap(:,:,itemp))
nbin = int(delta_max/(0.005d0/1000.d0))
delta_max = 1.1d0 * maxval(Agap(:,:,itemp))
nbin = NINT(delta_max / eps5) + 1
dbin = delta_max / dble(nbin)
IF ( .not. ALLOCATED(delta_k_bin) ) ALLOCATE( delta_k_bin(nbin) )
delta_k_bin(:) = zero
......@@ -2119,11 +2117,9 @@
DO ik = 1, nkfs
DO ibnd = 1, nbndfs
IF ( abs( ekfs(ibnd,ik) - ef0 ) .lt. fsthick ) THEN
DO ibin = 1, nbin
sigma = 1.d0 * dbin
weight = w0gauss( ( Agap(ibnd,ik,itemp) - dble(ibin) * dbin) / sigma, 0 ) / sigma
delta_k_bin(ibin) = delta_k_bin(ibin) + weight
ENDDO
ibin = nint( Agap(ibnd,ik,itemp) / dbin ) + 1
weight = w0g(ibnd,ik)
delta_k_bin(ibin) = delta_k_bin(ibin) + weight
ENDIF
ENDDO
ENDDO
......
!
! Copyright (C) 2010-2016 Samuel Ponce', Roxana Margine, Carla Verdi, Feliciano Giustino
! Copyright (C) 2016-2018 Samuel Ponce'
! Copyright (C) 2016-2019 Samuel Ponce', Roxana Margine, Feliciano Giustino
!
! This file is distributed under the terms of the GNU General Public
! License. See the file `LICENSE' in the root directory of the
! present distribution, or http://www.gnu.org/copyleft.gpl.txt .
!
!
!-----------------------------------------------------------------------
SUBROUTINE print_ibte( iqq, iq, totq, ef0, efcb, first_cycle, ind_tot, ind_totcb, &
......
This diff is collapsed.
......@@ -329,13 +329,13 @@
wgkq = wgauss( -ekq*inv_eptemp0, -99)
!
! = k-point weight * [f(E_k) - f(E_k+q)]/ [E_k+q - E_k -w_q + id]
! This is the imaginary part of the phonon self-energy, sans
! This is the imaginary part of minus the phonon self-energy, sans
! the matrix elements
!
!weight = wkf (ikk) * (wgkk - wgkq) * &
! aimag ( cone / ( ekq - ekk - wq - ci * degaussw0 ) )
!
! SP: The expression below (phonon self-energy) corresponds to
! SP: The expression below (minus phonon self-energy) corresponds to
! = pi*k-point weight*[f(E_k) - f(E_k+q)]*delta[E_k+q - E_k - w_q]
weight = pi * wkf (ikk) * (wgkk - wgkq)* &
w0gauss ( (ekq - ekk - wq) / degaussw0, 0) / degaussw0
......
This diff is collapsed.
......@@ -1423,7 +1423,7 @@
degaussw, nkf1, nkf2, nkf3
USE eliashbergcom, ONLY : nkfs, nbndfs, g2, ixkqf, ixqfs, nqfs, w0g, ekfs, ef0, dosef, wsph, &
wkfs, dwsph, a2f_iso, ixkff
USE constants_epw, ONLY : ryd2ev
USE constants_epw, ONLY : ryd2ev, eps2, zero, eps16
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : inter_pool_comm, my_pool_id, npool
USE mp_world, ONLY : mpime
......@@ -1599,23 +1599,25 @@
IF ( ALLOCATED(a2f) ) DEALLOCATE( a2f )
IF ( ALLOCATED(a2f_modeproj) ) DEALLOCATE( a2f_modeproj )
!
nbink = int( 1.25d0 * maxval(lambda_k(:,:)) / 0.005d0 )
dbink = 1.25d0 * maxval(lambda_k(:,:)) / dble(nbink)
nbink = NINT( 1.1d0 * MAXVAL(lambda_k(:,:)) / eps2 ) + 1
dbink = 1.1d0 * MAXVAL(lambda_k(:,:)) / DBLE(nbink)
!
IF ( .not. ALLOCATED(lambda_k_bin) ) ALLOCATE ( lambda_k_bin(nbink) )
lambda_k_bin(:) = 0.d0
lambda_k_bin(:) = zero
!
!SP : Should be initialized
nbin = 0
dbin = 0.0_DP
dbin = zero
!
IF ( iverbosity == 2 ) THEN
nbin = int( 1.25d0 * maxval(lambda_max(:)) / 0.005d0 )
dbin = 1.25d0 * maxval(lambda_max(:)) / dble(nbin)
IF ( .not. ALLOCATED(lambda_pairs) ) ALLOCATE ( lambda_pairs(nbin) )
lambda_pairs(:) = 0.d0
nbin = nint( 1.1d0 * MAXVAL(lambda_max(:)) / eps2 ) + 1
dbin = 1.1d0 * MAXVAL(lambda_max(:)) / dble(nbin)
IF ( .not. ALLOCATED(lambda_pairs) ) ALLOCATE ( lambda_pairs(nbin) )
lambda_pairs(:) = zero
ENDIF
!
WRITE(stdout,'(5x,a13,f21.7,a18,f21.7)') 'lambda_max = ', maxval(lambda_max(:)), ' lambda_k_max = ', maxval(lambda_k(:,:))
WRITE(stdout,'(5x,a13,f21.7,a18,f21.7)') 'lambda_max = ', maxval(lambda_max(:)), &
' lambda_k_max = ', maxval(lambda_k(:,:))
WRITE(stdout,'(a)') ' '
!
lambda_k(:,:) = 0.d0
......@@ -1631,20 +1633,16 @@
CALL lambdar_aniso_ver1( ik, iq, ibnd, jbnd, 0.d0, lambda_eph )
lambda_k(ik,ibnd) = lambda_k(ik,ibnd) + weight * lambda_eph
IF ( iverbosity == 2 ) THEN
DO ibin = 1, nbin
sigma = 1.d0 * dbin
weight = w0gauss( ( lambda_eph - dble(ibin) * dbin ) / sigma, 0 ) / sigma
lambda_pairs(ibin) = lambda_pairs(ibin) + weight
ENDDO
ibin = NINT( lambda_eph / dbin ) + 1
weight = w0g(ibnd,ik) * w0g(jbnd,ixkqf(ik,iq0))
lambda_pairs(ibin) = lambda_pairs(ibin) + weight
ENDIF
ENDIF
ENDDO ! jbnd
ENDDO ! iq
DO ibin = 1, nbink
sigma = 1.d0 * dbink
weight = w0gauss( ( lambda_k(ik,ibnd) - dble(ibin) * dbink ) / sigma, 0 ) / sigma
lambda_k_bin(ibin) = lambda_k_bin(ibin) + weight
ENDDO
ibin = NINT( lambda_k(ik,ibnd) / dbink ) + 1
weight = w0g(ibnd,ik)
lambda_k_bin(ibin) = lambda_k_bin(ibin) + weight
ENDIF
ENDDO ! ibnd
ENDDO ! ik
......@@ -1676,7 +1674,7 @@
OPEN(unit = iufillambda, file = TRIM(prefix)//".lambda_k_pairs", form = 'formatted')
WRITE(iufillambda,'(a12,a30)') '# lambda_nk',' \rho(lambda_nk) scaled to 1'
DO ibin = 1, nbink
WRITE(iufillambda,'(2f21.7)') dbink*dble(ibin), lambda_k_bin(ibin)/maxval(lambda_k_bin(:))
WRITE(iufillambda,'(2f21.7)') dbink*dble(ibin), lambda_k_bin(ibin)/MAXVAL(lambda_k_bin(:))
ENDDO
CLOSE(iufillambda)
!
......@@ -1761,9 +1759,4 @@
!
END SUBROUTINE evaluate_a2f_lambda
!
!
END MODULE superconductivity_aniso
!
! Copyright (C) 2010-2016 Samuel Ponce', Roxana Margine, Carla Verdi, Feliciano Giustino
! Copyright (C) 2016-2018 Samuel Ponce'
! Copyright (C) 2016-2019 Samuel Ponce', Roxana Margine, Feliciano Giustino
!
! This file is distributed under the terms of the GNU General Public
! License. See the file `LICENSE' in the root directory of the
......
......@@ -163,12 +163,17 @@
!! Used for the averaging
REAL(kind=DP) :: ekk2
!! Use for averaging
!
REAL(kind=DP) :: xkf_tmp (3, nkqtotf)
!! Temporary k-point coordinate (dummy variable)
REAL(kind=DP) :: wkf_tmp(nkqtotf)
!! Temporary k-weights (dummy variable)
REAL(kind=DP) :: dfnk
!! df/de
REAL(kind=DP) :: etemp
!! Temperature