1! ---
2! Copyright (C) 1996-2016	The SIESTA group
3!  This file is distributed under the terms of the
4!  GNU General Public License: see COPYING in the top directory
5!  or http://www.gnu.org/copyleft/gpl.txt .
6! See Docs/Contributors.txt for a list of contributors.
7! ---
8      module m_projected_DOS
9
10      use precision
11
12      implicit none
13
14      private
15
16      public :: init_projected_DOS, projected_DOS
17
18      logical  :: different_pdos_grid    ! Indicates if the grid is the same as the SCF one or not
19
20      CONTAINS
21
22      subroutine init_projected_DOS()
23
24      USE siesta_options
25      use fdf,         only: fdf_block, fdf_isblock
26      use Kpoint_pdos
27      use parallel,    only: IOnode, Nodes
28      use siesta_geom, only: ucell
29      use m_spin,      only: spin
30
31!-------------------------------------------------------------------------BEGIN
32!     Compute the projected density of states
33      do_pdos = fdf_isblock('ProjectedDensityOfStates')
34      if ( .not. do_pdos ) return
35
36      if (isolve.ne.SOLVE_DIAGON) then
37         if (.not.((isolve.eq.SOLVE_MINIM).and.
38     .        minim_calc_eigenvalues)) then
39            if (IONode) then
40               write(6,*)
41     .              'siesta: ERROR: PDOS implemented only with diagon'
42            endif
43            do_pdos = .false.
44         endif
45      endif
46
47      if ( .not. do_pdos ) return
48
49!     Determine whether the projected density of states is to be computed
50!     on a different grid to the SCF calculation
51      call setup_Kpoint_pdos( ucell, different_pdos_grid )
52
53!---------------------------------------------------------------------------END
54
55      end subroutine init_projected_DOS
56
57      subroutine projected_DOS()
58
59      use sparse_matrices
60      USE siesta_options
61      use alloc,       only : re_alloc
62      use atomlist,    only : indxuo, no_s, no_u, no_l
63      use fdf
64      use sys,         only : die
65      use Kpoint_grid
66      use Kpoint_pdos
67      use parallel,    only: IOnode
68      use m_eo
69      use m_spin,      only: h_spin_dim, spinor_dim
70      use units,       only: eV
71
72      implicit none
73
74      type(block_fdf)            :: bfdf
75      type(parsed_line), pointer :: pline
76
77      real(dp) :: factor
78      logical  :: dummy ! Logical to hold return value from call to fdf_block
79      integer  :: nhist ! Number of histogram intervals in projected DOS
80      real(dp) :: e1    ! Lower bound of energy range
81      real(dp) :: e2    ! Upper bound of energy range
82      real(dp) :: sigma ! Energy width used to convolute partial DOS
83
84!------------------------------------------------------------------------- BEGIN
85! Compute the projected density of states
86
87      if ( .not. do_PDOS ) return
88
89      ! Call fdf_block to get iu - presence has already been tested in init_projected_PDOS
90      dummy = fdf_block('ProjectedDensityOfStates',bfdf)
91! Find the desired energy range
92      if (.not. fdf_bline(bfdf,pline))
93     $     call die('projected_DOS: ERROR in ' //
94     $     'ProjectedDensityOfStates block')
95      if ((fdf_bnreals(pline) .lt. 3) .or.
96     $     (fdf_bnnames(pline) .ne. 1))
97     $   call die("Wrong format in PDOS block, not enough reals/names")
98      factor = fdf_convfac( fdf_bnames(pline,1), 'Ry' )
99      e1 = fdf_breals(pline,1) * factor
100      e2 = fdf_breals(pline,2) * factor
101      sigma = fdf_breals(pline,3) * factor
102      nhist = fdf_bintegers(pline,1)
103      if (IOnode) then
104         write(6,'(a)') 'siesta: PDOS info: '
105         write(6,'(a,3(f8.2,a),2x,i5)')
106     $        'siesta: e1, e2, sigma, nhist: ',
107     $        e1/eV,' eV',e2/eV,' eV',sigma/eV,' eV', nhist
108      end if
109      call fdf_bclose(bfdf)
110
111! If the k points have been set specifically for the PDOS then use this set
112      if (different_pdos_grid) then
113
114! If the number of k points has increased then reallocate eo and qo
115         if (maxk_pdos.gt.nkpnt) then
116            call re_alloc(eo,1,no_u,1,spinor_dim,1,maxk_pdos,name="eo",
117     .           routine="projected_dos")
118         endif
119
120         call pdos( no_s, h_spin_dim, spinor_dim, no_l,
121     .        maxnh,
122     .        no_u, numh, listhptr, listh, H, S,
123     .        e1, e2, sigma, nhist, xijo, indxuo, gamma_PDOS,
124     .        nkpnt_pdos, kpoints_pdos, kweight_pdos, eo,
125     .        no_u)
126      else
127! otherwise use the SCF grid
128         call pdos( no_s, h_spin_dim, spinor_dim, no_l,
129     .        maxnh,
130     .        no_u, numh, listhptr, listh, H, S,
131     .        e1, e2, sigma, nhist,
132     .        xijo, indxuo, gamma_SCF,
133     .        nkpnt, kpoint, kweight, eo,
134     .        no_u)
135      endif
136
137      end subroutine projected_DOS
138
139      end module m_projected_DOS
140