Hello developers,
There is a bug in the way SIESTA reads the ProjectedDensityOfStates
block
which can cause the number of energy points to be read
incorrectly. The
energy range and broadening are read as the first three values on
the line,
but the number of energy points is read as the first integer on the
line.
Therefore, if any other value of the line is written as an integer,
then
that will be incorrectly taken as the number of energy points. If
this
value happens to be negative a number of MPI errors are generated
as the
code attempts to deal with an array with a negative size, which is
how I
noticed this problem. The code should be changed to read the
fourth value
on the line, not the first integer, so that it conforms to the
specifications in the manual.
Consider below the following pairs of FDF input and standard output
showing
what I mean. The first example is correct behavior; the rest are
wrong.
Also at the bottom is the section of code that reads and writes the
values
from the PDOS block, and I believe "nhist = integers(p,1)" is the
line that
needs to be fixed.
David Strubbe
UC Berkeley
===============
%block ProjectedDensityOfStates
20.0 30.0 0.01 10 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: 20.00 eV 30.00 eV 0.01 eV
10
---
%block ProjectedDensityOfStates
-30 -0 0.01 6000 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: -30.00 eV 0.00 eV 0.01 eV
-30
---
%block ProjectedDensityOfStates
-30 -20 0.01 100 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: -30.00 eV -20.00 eV 0.01 eV
-30
---
%block ProjectedDensityOfStates
-20 10 0.01 100 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: -20.00 eV 10.00 eV 0.01 eV
-20
---
%block ProjectedDensityOfStates
20 30 0.05 10 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: 20.00 eV 30.00 eV 0.05 eV
20
---
%block ProjectedDensityOfStates
20 30 0.01 10 eV
%endblock ProjectedDensityOfStates
siesta: PDOS info:
siesta: e1, e2, sigma, nhist: 20.00 eV 30.00 eV 0.01 eV
20
============================
siesta.F:
C Find the desired energy range
if (IOnode) then
read(iu,'(a)') line
p=>digest(line)
if (nvalues(p).lt.3 .or. nnames(p).ne.1)
$ call die("Wrong format in PDOS block")
factor = fdf_convfac( names(p,1), 'Ry' )
e1 = values(p,1) * factor
e2 = values(p,2) * factor
sigma = values(p,3) * factor
nhist = integers(p,1)
write(6,'(a)') 'siesta: PDOS info: '
write(6,'(a,3(f8.2,a),2x,i5)')
$ 'siesta: e1, e2, sigma, nhist: ',
$ e1/eV,' eV',e2/eV,' eV',sigma/eV,' eV', nhist
endif
