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
