Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 : MODULE m_types_eigdos
7 : USE m_juDFT
8 : use m_constants
9 : IMPLICIT NONE
10 : PRIVATE
11 : PUBLIC:: t_eigdos,t_eigdos_list,t_eigdos_make_dos
12 :
13 : TYPE,abstract :: t_eigdos
14 : CHARACTER(len=20) :: name_of_dos="unnamed"
15 : !each eigenvalue might be described by weights
16 : REAL,ALLOCATABLE :: eig(:,:,:)
17 : REAL,ALLOCATABLE :: dos_grid(:) !This is the grid the DOS part uses internally (FOR IO use the routine below)
18 : REAL,ALLOCATABLE :: dos(:,:,:) !(grid,spins,weights)
19 : CONTAINS
20 : procedure(get_weight_name),DEFERRED :: get_weight_name
21 : procedure(get_num_weights),DEFERRED :: get_num_weights
22 : procedure(get_weight_eig),DEFERRED :: get_weight_eig !should be overwritten in derived type
23 : procedure :: get_spins
24 : procedure :: get_neig
25 : procedure :: get_eig
26 : procedure :: get_dos_grid
27 : procedure :: make_dos=>t_eigdos_make_dos
28 : procedure :: smooth=>dosdata_smooth
29 : procedure :: write_raw !should be implemented later to allow eig66 functionality
30 : procedure :: write_dos
31 : procedure :: write_band
32 : procedure :: write_EVData
33 : procedure :: sym_weights,sym_weights_eigdos
34 : END TYPE
35 :
36 : type::t_eigdos_list
37 : class(t_eigdos),POINTER:: p
38 : end type
39 :
40 : INTERFACE
41 : function get_weight_eig(this,id)
42 : import t_eigdos
43 : class(t_eigdos),intent(in):: this
44 : INTEGER,intent(in) :: id
45 : real,allocatable:: get_weight_eig(:,:,:)
46 : end function
47 : END interface
48 : INTERFACE
49 : integer function get_num_weights(this)
50 : import t_eigdos
51 : class(t_eigdos),intent(in):: this
52 : end function
53 : end interface
54 : INTERFACE
55 : character(len=20) function get_weight_name(this,id)
56 : import t_eigdos
57 : class(t_eigdos),intent(in):: this
58 : INTEGER,intent(in) :: id
59 : end function
60 : end interface
61 : CONTAINS
62 1 : subroutine sym_weights(this)
63 : class(t_eigdos),INTENT(INOUT):: this
64 1 : end subroutine
65 :
66 142 : subroutine sym_weights_eigdos(this,weight)
67 : !Subroutine to make all weights equal for degenerate states
68 : class(t_eigdos),INTENT(IN):: this
69 : real,intent(inout) :: weight(:,:,:)
70 :
71 : INTEGER:: ispin,ikpt,i,j
72 :
73 304 : DO ispin=1,size(this%eig,3)
74 3092 : DO ikpt=1,size(this%eig,2)
75 : ! Make sure equivalent states have same weight
76 : i=1
77 53868 : DO while(i<size(this%eig,1))
78 : j=1
79 68300 : do while (abs(this%eig(i,ikpt,ispin)-this%eig(i+j,ikpt,ispin))<1E-5)
80 18850 : j=j+1
81 68300 : if (i+j>size(this%eig,1)) exit
82 : ENDDO
83 50918 : if (j>1) THEN
84 16864 : j=j-1
85 : !Make sure all weights are equal
86 88292 : weight(i:i+j,ikpt,ispin)=sum(weight(i:i+j,ikpt,ispin))/(j+1)
87 : i=i+j
88 : endif
89 50918 : i=i+1
90 : enddo
91 : ENDDO
92 : ENDDO
93 142 : end subroutine
94 :
95 71 : function get_neig(this)
96 : CLASS(t_eigdos),INTENT(IN)::this
97 : integer,allocatable :: get_neig(:,:)
98 71 : real,allocatable::ev(:,:,:)
99 : integer ::k,j
100 :
101 30252 : ev=this%get_eig()
102 284 : allocate(get_neig(size(ev,2),size(ev,3)))
103 144 : DO j=1,this%get_spins()
104 1572 : DO k=1,size(ev,2)
105 30181 : get_neig(k,j)=count(ev(:,k,j)<1E99)
106 : ENDDO
107 : ENDDO
108 : end function
109 :
110 158 : pure integer function get_spins(this)
111 : CLASS(t_eigdos),INTENT(IN)::this
112 158 : get_spins=size(this%eig,3)
113 158 : END function
114 528 : function get_eig(this)
115 : CLASS(t_eigdos),INTENT(IN):: this
116 : real,allocatable:: get_eig(:,:,:)
117 123219 : get_eig=this%eig
118 : END function
119 :
120 184 : function get_dos_grid(this)
121 : CLASS(t_eigdos),INTENT(IN):: this
122 : real,allocatable:: get_dos_grid(:)
123 121900 : get_dos_grid=this%dos_grid
124 : END function
125 :
126 5 : subroutine dosdata_smooth(eigdos,banddos)
127 : use m_smooth
128 : use m_types_banddos
129 : class(t_eigdos),INTENT(INOUT) :: eigdos
130 : type(t_banddos),INTENT(IN) :: banddos
131 :
132 : integer :: jspin,i
133 5 : real,allocatable :: dos_grid(:)
134 :
135 5 : if (.not.allocated(eigdos%dos)) return
136 20 : if (size(eigdos%dos)==0) return
137 20 : if (size(eigdos%dos)<1) return
138 5 : if (banddos%sig_dos==0.0) RETURN
139 6610 : dos_grid=eigdos%get_dos_grid()
140 :
141 13220 : IF ( NINT((maxval(dos_grid) - minval(dos_grid))/banddos%sig_dos) > size(dos_grid) ) THEN
142 0 : WRITE(oUnit,*) 'sig_dos too small for DOS smoothing:'
143 0 : WRITE(oUnit,*) 'Reduce energy window or enlarge banddos%sig_dos!'
144 0 : WRITE(oUnit,*) 'For now: no smoothing done'
145 0 : return
146 : ENDIF
147 :
148 87 : DO i=1,size(eigdos%dos,3)
149 176 : DO jspin=1,eigdos%get_spins()
150 171 : CALL smooth(dos_grid,eigdos%dos(:,jspin,i),banddos%sig_dos,size(eigdos%dos_grid))
151 : ENDDO
152 : ENDDO
153 5 : END subroutine
154 :
155 5 : subroutine write_dos(eigdos,hdf_id)
156 : #ifdef CPP_HDF
157 : use HDF5
158 : use m_banddos_io
159 : #endif
160 : class(t_eigdos),INTENT(INOUT):: eigdos
161 : #ifdef CPP_HDF
162 : integer(HID_T),intent(in) ::hdf_id
163 : #else
164 : integer, intent(in) ::hdf_id
165 : #endif
166 : integer:: jspin,i,ind,id, n
167 : character(len=100)::filename
168 5 : real,allocatable:: dos_grid(:)
169 : LOGICAL l_printTextDOS
170 :
171 5 : l_printTextDOS = .TRUE.
172 :
173 : #ifdef CPP_HDF
174 87 : DO n=1,eigdos%get_num_weights()
175 82 : print *, "writedos:",n,eigdos%get_num_weights()
176 87 : call writedosData(hdf_ID,eigdos%name_of_dos,eigdos%get_dos_grid(),eigdos%get_weight_name(n),eigdos%dos(:,:,n))
177 : enddo
178 5 : IF(eigdos%get_num_weights().GT.40) THEN
179 1 : WRITE(*,*) 'Number of weights in ', TRIM(ADJUSTL(eigdos%name_of_dos)),' DOS too large for simple text output.'
180 1 : WRITE(*,*) 'Output only in banddos.hdf file.'
181 : l_printTextDOS = .FALSE.
182 : END IF
183 : #endif
184 :
185 1 : IF (.NOT.l_printTextDOS) RETURN
186 4 : if (.not.allocated(eigdos%dos)) return
187 16 : if (size(eigdos%dos)==0) return
188 9 : DO jspin=1,eigdos%get_spins()
189 5 : write(filename,"(a,a,i0)") trim(eigdos%name_of_dos),".",jspin
190 5 : open(999,file=filename)
191 48 : write(999,"(999a21)") "#energy",(eigdos%get_weight_name(id),id=1,eigdos%get_num_weights())
192 48 : write(*,"(999a21)") filename,(eigdos%get_weight_name(id),id=1,eigdos%get_num_weights())
193 6615 : dos_grid=eigdos%get_dos_grid()
194 6610 : DO i=1,size(dos_grid)
195 63413 : write(999,"(999(e20.8,1x))") dos_grid(i)*hartree_to_ev_const,(eigdos%dos(i,jspin,id)/hartree_to_ev_const,id=1,eigdos%get_num_weights())
196 : ENDDO
197 5 : close(999)
198 9 : write(*,*) "done:",filename
199 : ENDDO
200 4 : END subroutine
201 :
202 2 : subroutine write_band(eigdos,kpts,title,cell,hdf_id,efermi,banddos)
203 : use m_types_kpts
204 : use m_types_cell
205 : use m_gnuplot_BS
206 : use m_types_banddos
207 : #ifdef CPP_HDF
208 : use HDF5
209 : use m_banddos_io
210 : #endif
211 : class(t_eigdos),INTENT(INOUT):: eigdos
212 : type(t_kpts),intent(in) :: kpts
213 : type(t_banddos),INTENT(IN) :: banddos
214 : type(t_cell),intent(in) :: cell
215 : CHARACTER(LEN=*), INTENT(IN) :: title
216 : real,intent(in) :: efermi
217 : #ifdef CPP_HDF
218 : integer(HID_T),intent(in) ::hdf_id
219 : INTEGER::n
220 :
221 : #else
222 : integer,intent(in):: hdf_id !not used
223 : #endif
224 :
225 : integer:: jspin,i,k
226 2 : real,allocatable :: ev(:,:,:),kx(:)
227 : real :: vkr(3),vkr_prev(3)
228 : character(len=100)::file
229 6 : allocate(kx(kpts%nkpt))
230 : #ifdef CPP_HDF
231 28 : DO n=1,eigdos%get_num_weights()
232 28 : call writebandData(hdf_id,kpts,eigdos%name_of_dos,eigdos%get_weight_name(n),eigdos%get_weight_eig(n),eigdos%get_eig())
233 : enddo
234 : #endif
235 2 : if (eigdos%name_of_dos.ne."Local") then
236 : #ifndef CPP_HDF
237 : print *,"WARNING, only very basic BS written without HDF"
238 : #endif
239 0 : return
240 : endif
241 5 : DO jspin=1,eigdos%get_spins()
242 3 : write(file,"(a,i0)") "bands.",jspin
243 3 : open(18,file=file)
244 4279 : ev=eigdos%get_eig()
245 3 : kx(1) = 0.0
246 39 : vkr_prev=matmul(kpts%bk(:,1),cell%bmat)
247 44 : DO k = 2, kpts%nkpt
248 533 : vkr=matmul(kpts%bk(:,k),cell%bmat)
249 164 : kx(k)=kx(k-1)+ sqrt(dot_product(vkr-vkr_prev,vkr-vkr_prev))
250 44 : vkr_prev=vkr
251 : ENDDO
252 254 : DO i = 1, minval(eigdos%get_neig())
253 2461 : DO k = 1, kpts%nkpt
254 2458 : write(18,'(999f15.9)') kx(k),(ev(i,k,jspin)-efermi)*hartree_to_ev_const
255 : !-eFermiCorrection
256 : ENDDO
257 : ENDDO
258 5 : CLOSE (18)
259 : enddo
260 2 : call gnuplot_bs(kpts,title,cell,eigdos%get_spins())
261 2 : IF (banddos%unfoldband) call write_gnu_sc(banddos,kpts,title,cell,eigdos%get_spins())
262 2 : end subroutine
263 :
264 5 : subroutine write_EVData(eigdos,hdf_id)
265 : #ifdef CPP_HDF
266 : use HDF5
267 : use m_banddos_io
268 : #endif
269 : class(t_eigdos),INTENT(INOUT):: eigdos
270 : #ifdef CPP_HDF
271 : integer(HID_T),intent(in) ::hdf_id
272 : INTEGER::n
273 : #else
274 : integer,intent(in):: hdf_id !not used
275 : #endif
276 :
277 : #ifdef CPP_HDF
278 87 : DO n = 1, eigdos%get_num_weights()
279 87 : CALL writeEVData(hdf_id,eigdos%name_of_dos,eigdos%get_weight_name(n),eigdos%get_eig(),eigdos%get_weight_eig(n))
280 : END DO
281 : #endif
282 5 : end subroutine
283 :
284 5 : subroutine t_eigdos_make_dos(eigdos,kpts,input,banddos,efermi)
285 : use m_types_banddos
286 : use m_types_input
287 : use m_dosbin
288 : use m_ptdos
289 : use m_tetra_dos
290 : use m_dostetra
291 : use m_types_kpts
292 :
293 : class(t_eigdos),intent(inout):: eigdos
294 : type(t_banddos),intent(in) :: banddos
295 : type(t_kpts),intent(in) :: kpts
296 : type(t_input),intent(in) :: input
297 : real,intent(in) :: efermi
298 :
299 : integer ::n
300 : real :: emin,emax
301 :
302 5 : emin=min(banddos%e1_dos,banddos%e2_dos)-efermi
303 5 : emax=max(banddos%e1_dos,banddos%e2_dos)-efermi
304 5 : if (allocated(eigdos%dos)) return
305 :
306 25 : allocate(eigdos%dos(banddos%ndos_points,eigdos%get_spins(),eigdos%get_num_weights()))
307 : !Generate DOS grid
308 5 : if (allocated(eigdos%dos_grid)) deallocate(eigdos%dos_grid)
309 15 : allocate(eigdos%dos_grid(banddos%ndos_points))
310 6610 : DO n=1,banddos%ndos_points
311 6610 : eigdos%dos_grid(n)=emin+(emax-emin)/(banddos%ndos_points-1.0)*(n-1.0)
312 : ENDDO
313 :
314 87 : DO n=1,eigdos%get_num_weights()
315 82 : print *,eigdos%name_of_dos,n,eigdos%get_num_weights()
316 5 : SELECT CASE(input%bz_integration)
317 :
318 : CASE(BZINT_METHOD_HIST, BZINT_METHOD_GAUSS)
319 14 : call dos_bin(input%jspins,kpts%wtkpt,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
320 : CASE(BZINT_METHOD_TRIA)
321 68 : IF(input%film) THEN
322 0 : CALL ptdos(input%jspins,kpts,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
323 : ELSE
324 : CALL tetra_dos(input%jspins,kpts,eigdos%dos_grid,eigdos%get_neig(),eigdos%get_eig(),&
325 68 : eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
326 : ENDIF
327 : CASE(BZINT_METHOD_TETRA)
328 82 : CALL dostetra(kpts,input,eigdos%dos_grid,eigdos%get_eig(),eigdos%get_weight_eig(n),eigdos%dos(:,:,n),efermi)
329 : END SELECT
330 :
331 : end do
332 : END subroutine
333 :
334 :
335 :
336 0 : subroutine write_raw(this,id)
337 : class(t_eigdos),INTENT(IN):: this
338 : INTEGER,INTENT(IN) :: id
339 :
340 :
341 0 : end subroutine
342 :
343 :
344 :
345 0 : END MODULE m_types_eigdos
|