unfoldBandKPTS.f90 22.9 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7 8 9 10
MODULE m_unfold_band_kpts

CONTAINS

11
  SUBROUTINE build_primitive_cell(banddos,p_cell,cell)
12 13 14 15 16 17 18
    USE m_types
    USE m_inv3
    USE m_constants, ONLY : tpi_const
    implicit none
    TYPE(t_banddos),INTENT(IN)  :: banddos
    TYPE(t_cell),INTENT(IN)     :: cell
    TYPE(t_cell),INTENT(INOUT)  :: p_cell
19

20 21 22 23 24
    INTEGER :: i
    DO i =1,3
	p_cell%amat(1,i)=cell%amat(1,i)/banddos%s_cell_x
	p_cell%amat(2,i)=cell%amat(2,i)/banddos%s_cell_y 
	p_cell%amat(3,i)=cell%amat(3,i)/banddos%s_cell_z 
25 26 27
!	p_cell%amat(i,1)=cell%amat(i,1)/banddos%s_cell_x
!	p_cell%amat(i,2)=cell%amat(i,2)/banddos%s_cell_y 
!	p_cell%amat(i,3)=cell%amat(i,3)/banddos%s_cell_z 
28 29 30 31
    END DO
    CALL inv3(p_cell%amat,p_cell%bmat,p_cell%omtil)
    p_cell%bmat=p_cell%bmat*tpi_const
    p_cell%latnam=cell%latnam
32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47
  END SUBROUTINE  build_primitive_cell

  SUBROUTINE unfold_band_kpts(banddos,p_cell,cell,p_kpts,kpts)
    USE m_types
    USE m_inv3
    USE m_constants, ONLY : tpi_const

    implicit none

    TYPE(t_banddos),INTENT(IN)  :: banddos
    TYPE(t_cell),INTENT(IN)     :: cell
    TYPE(t_cell),INTENT(INOUT)  :: p_cell
    TYPE(t_kpts),INTENT(INOUT)  :: p_kpts
    TYPE(t_kpts),INTENT(INOUT)  :: kpts
   
    CALL build_primitive_cell(banddos,p_cell,cell)
48 49

    p_kpts=kpts
50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69
    !write(1088,*) 'banddos%unfoldband: ', banddos%unfoldband
    !write(1088,*) 'brav. matrix: '
    !write(1088,'(f15.8,f15.8,f15.8)') cell%amat(1,1), cell%amat(1,2), cell%amat(1,3)
    !write(1088,'(f15.8,f15.8,f15.8)') cell%amat(2,1), cell%amat(2,2), cell%amat(2,3)
    !write(1088,'(f15.8,f15.8,f15.8)') cell%amat(3,1), cell%amat(3,2), cell%amat(3,3)
    !write(1088,*) 'brav. rez. matrix: '
    !write(1088,'(f15.8,f15.8,f15.8)') cell%bmat(1,1), cell%bmat(1,2), cell%bmat(1,3)
    !write(1088,'(f15.8,f15.8,f15.8)') cell%bmat(2,1), cell%bmat(2,2), cell%bmat(2,3)
    !write(1088,'(f15.8,f15.8,f15.8)') cell%bmat(3,1), cell%bmat(3,2), cell%bmat(3,3)
    !write(1088,*) ' primitive brav. matrix: '
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%amat(1,1), p_cell%amat(1,2), p_cell%amat(1,3)
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%amat(2,1), p_cell%amat(2,2), p_cell%amat(2,3)
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%amat(3,1), p_cell%amat(3,2), p_cell%amat(3,3)
    !write(1088,*) 'primitive brav. rez. matrix: '
    !write(89,'(3f15.8)') p_cell%bmat
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%bmat(1,1), p_cell%bmat(1,2), p_cell%bmat(1,3)
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%bmat(2,1), p_cell%bmat(2,2), p_cell%bmat(2,3)
    !write(1088,'(f15.8,f15.8,f15.8)') p_cell%bmat(3,1), p_cell%bmat(3,2), p_cell%bmat(3,3)
    !write(1088,'(a,i7,a,i7)') 'kpts%nkpt',kpts%nkpt,'   p_kpts%nkpt',p_kpts%nkpt
    !write(1088,*) kpts%specialPoints
70 71 72 73 74
  END SUBROUTINE unfold_band_kpts
  
  SUBROUTINE find_supercell_kpts(banddos,p_cell,cell,p_kpts,kpts)
    USE m_types
    USE m_juDFT
75
    USE m_inv3
76 77 78 79
    implicit none

    TYPE(t_banddos),INTENT(IN)  :: banddos
    TYPE(t_cell),INTENT(IN)     :: cell
80 81
    TYPE(t_cell),INTENT(IN)     :: p_cell
    TYPE(t_kpts),INTENT(IN)     :: p_kpts
82
    TYPE(t_kpts),INTENT(INOUT)  :: kpts
83 84
    
    INTEGER :: i,m1,m2,m3
85 86
    REAL    :: rez_inv_to_internal(3,3)
    REAL    :: rez_inv_det
87
    REAL    :: list(13,p_kpts%nkpt)  !cartesion coordinates for k,K,m
88 89 90 91
    REAL    :: pc_kpoint_i(3)    !primitive cell kpoint internal
    REAL    :: sc_kpoint_i(3)    !super cell kpoint internal
    REAL    :: pc_kpoint_c(3)    !primitive cell kpoint cartesian
    REAL    :: sc_kpoint_c(3)    !super cell kpoint cartesian
92
    REAL    :: eps(3)
93
    REAL    :: eps_r, eps_kpt
94
    LOGICAL :: representation_found
95
    REAL    ::kpt_dist
96

97
    eps = 1.0e-10
98
    eps_r = 0.000000001
99 100

    CALL inv3(cell%bmat,rez_inv_to_internal,rez_inv_det)
101 102
    !write(1088,*) p_kpts%specialPoints
    !write(333,'(3f15.8)')p_kpts%bk
103
    kpt_dist=0
104
    DO i= 1,size(list,2)
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150
	!        pc_kpoint_c(1)=p_kpts%bk(1,i)*p_cell%bmat(1,1)+p_kpts%bk(2,i)*p_cell%bmat(1,2)+p_kpts%bk(3,i)*p_cell%bmat(1,3)
	!        pc_kpoint_c(2)=p_kpts%bk(1,i)*p_cell%bmat(2,1)+p_kpts%bk(2,i)*p_cell%bmat(2,2)+p_kpts%bk(3,i)*p_cell%bmat(2,3)
	!        pc_kpoint_c(3)=p_kpts%bk(1,i)*p_cell%bmat(3,1)+p_kpts%bk(2,i)*p_cell%bmat(3,2)+p_kpts%bk(3,i)*p_cell%bmat(3,3)
		pc_kpoint_c(1)=p_kpts%bk(1,i)*p_cell%bmat(1,1)+p_kpts%bk(2,i)*p_cell%bmat(2,1)+p_kpts%bk(3,i)*p_cell%bmat(3,1)
		pc_kpoint_c(2)=p_kpts%bk(1,i)*p_cell%bmat(1,2)+p_kpts%bk(2,i)*p_cell%bmat(2,2)+p_kpts%bk(3,i)*p_cell%bmat(3,2)
		pc_kpoint_c(3)=p_kpts%bk(1,i)*p_cell%bmat(1,3)+p_kpts%bk(2,i)*p_cell%bmat(2,3)+p_kpts%bk(3,i)*p_cell%bmat(3,3)
		list(1,i)=pc_kpoint_c(1)
		list(2,i)=pc_kpoint_c(2)
		list(3,i)=pc_kpoint_c(3)
	!!!!------- finding kpts in primitive rez. unit cell ----- 
	!	representation_found=.false.
	!m_loop:	DO m1= -banddos%s_cell_x,banddos%s_cell_x
	!		DO m2= -banddos%s_cell_y,banddos%s_cell_y
	!			DO m3= -banddos%s_cell_z,banddos%s_cell_z
	!				pc_kpoint_c(1)=list(1,i)-m1*cell%bmat(1,1)-m2*cell%bmat(1,2)-m3*cell%bmat(1,3)
	!				pc_kpoint_c(2)=list(2,i)-m1*cell%bmat(2,1)-m2*cell%bmat(2,2)-m3*cell%bmat(2,3)
	!				pc_kpoint_c(3)=list(3,i)-m1*cell%bmat(3,1)-m2*cell%bmat(3,2)-m3*cell%bmat(3,3)
	!!				IF (         (dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,1)) >= 0).AND.((dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,1)) < dot_product(cell%bmat(:,1), cell%bmat(:,1)))) &
	!!				     & .AND. (dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,2)) >= 0).AND.((dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,2)) < dot_product(cell%bmat(:,2), cell%bmat(:,2)))) &
	!!				     & .AND. (dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,3)) >= 0).AND.((dot_product(pc_kpoint_c(:)+eps(:), cell%bmat(:,3)) < dot_product(cell%bmat(:,3), cell%bmat(:,3))))) THEN
	!				IF (all((matmul(rez_inv_to_internal,pc_kpoint_c)+eps(:))>=0).and.all((matmul(rez_inv_to_internal,pc_kpoint_c)+eps(:))<1)) THEN
	!					list(4,i)=pc_kpoint_c(1)
	!					list(5,i)=pc_kpoint_c(2)
	!					list(6,i)=pc_kpoint_c(3)
	!					list(7,i)=-m1
	!					list(8,i)=-m2
	!					list(9,i)=-m3
	!					representation_found=.true.
	!				END IF
	 !      			        IF (representation_found) EXIT m_loop
	!			END DO
	!		END DO
	!	END DO m_loop
	 !       IF (.not.representation_found) THEN
	  !      write(*,'(a,f15.8,f15.8,f15.8)') 'No representation found for the following kpoint:',list(1,i),list(2,i),list(3,i)
	   !     END IF
	   !----------------------- method internal coordintes --------------------
	    sc_kpoint_i(:)=matmul(pc_kpoint_c,rez_inv_to_internal)
	    pc_kpoint_i(:)=p_kpts%bk(1:3,i)
	    !sc_kpoint_i(:) = sc_kpoint_i(:) + 0.5
	    m1 = FLOOR(sc_kpoint_i(1))
	    m2 = FLOOR(sc_kpoint_i(2))
	    m3 = FLOOR(sc_kpoint_i(3))
	    m1=0
	    m2=0
	    m3=0
151 152 153
	    sc_kpoint_i(1) = sc_kpoint_i(1) - m1
	    sc_kpoint_i(2) = sc_kpoint_i(2) - m2
	    sc_kpoint_i(3) = sc_kpoint_i(3) - m3 
154 155 156 157 158 159 160
	    !sc_kpoint_i(:) = sc_kpoint_i(:) - 0.5
	    list(4,i)=sc_kpoint_i(1)
	    list(5,i)=sc_kpoint_i(2)
	    list(6,i)=sc_kpoint_i(3)
	    list(7,i)=m1
	    list(8,i)=m2
	    list(9,i)=m3 !this whole block is to move kpoints into first BZ within -0.5 to 0.5
161

162
	!  	kpts%bk(:,i)=matmul(rez_inv_to_internal,pc_kpoint_c)
163 164 165
	    !-------------saving old kpts----------
	    list(11:13,i)=kpts%bk(:,i)
  	    !------finished---------
166 167 168 169 170 171
	    kpts%bk(:,i)=list(4:6,i)
	
	IF (i>1) THEN
	kpt_dist=kpt_dist+sqrt(dot_product(list(1:3,i)-list(1:3,i-1),list(1:3,i)-list(1:3,i-1)))
	END IF
	list(10,i)=kpt_dist
172
    END DO
173 174
    !write(91,'(3f15.8)') kpts%bk
    !write(92,*) kpts%wtkpt
175
    ALLOCATE (kpts%sc_list(13,p_kpts%nkpt))
176
    kpts%specialPointIndices(:) = p_kpts%specialPointIndices(:)
177
    kpts%sc_list=list
178
    !write(90,'(10f15.8)') kpts%sc_list
179 180
  END SUBROUTINE find_supercell_kpts

181
 SUBROUTINE calculate_plot_w_n(banddos,cell,kpts,smat_unfold,zMat,lapw,i_kpt,jsp,eig,results,input,atoms,unfoldingBuffer,mpi)
182 183 184 185 186 187 188
	USE m_types
	USE m_juDFT
	USE m_inv3
	USE m_types_mpimat
        USE m_constants
	implicit none

189 190
        TYPE(t_input),INTENT(IN) :: input
        TYPE(t_atoms),INTENT(IN)     :: atoms
191
	TYPE(t_banddos),INTENT(IN)  :: banddos
192
	TYPE(t_results),INTENT(INOUT)  :: results
193
	TYPE(t_cell),INTENT(IN)     :: cell
194
	TYPE(t_kpts),INTENT(INOUT)     :: kpts
195 196 197
	CLASS(t_mat),INTENT(INOUT)  :: smat_unfold
	CLASS(t_mat),INTENT(IN)     :: zMat
	TYPE(t_lapw),INTENT(IN)     :: lapw
198
        TYPE(t_mpi),INTENT(IN)       :: mpi
199 200 201
	TYPE(t_cell)      :: p_cell
	INTEGER, INTENT(IN)	    :: i_kpt,jsp
	REAL, INTENT(IN)	    :: eig(:)
202
        COMPLEX, INTENT(INOUT)         :: unfoldingBuffer(:,:,:)
203
	INTEGER :: i,j,k,l,n
204
	INTEGER :: na,n_i,nn,nk,nki,gi,lo
205
	REAL, ALLOCATABLE	::w_n(:)
206 207 208 209
	COMPLEX, ALLOCATABLE    ::w_n_c(:)
	REAL, ALLOCATABLE	::w_n_sum(:)
	COMPLEX, ALLOCATABLE    ::w_n_c_sum(:)
        LOGICAL :: method_rubel=.false.
210 211
        LOGICAL :: write_to_file=.false.
        CLASS(t_mat), ALLOCATABLE :: zMat_s
212

213
!	method_rubel=.true.    !this switch is to switch between overlap matrix and rubel method (without overlap matrix)
214

215 216 217 218 219 220 221 222 223 224 225 226
	CALL build_primitive_cell(banddos,p_cell,cell)
	IF (.not. method_rubel) THEN
		DO j = 1, lapw%nv(jsp)
		  DO i = 1, j-1
	      		IF(smat_unfold%l_real) THEN
				smat_unfold%data_r(j,i) = smat_unfold%data_r(i,j)
	      		ELSE
				smat_unfold%data_c(j,i) = CONJG(smat_unfold%data_c(i,j))
	      		END IF
		   END DO
		END DO
	END IF
227 228 229 230 231 232
!   	write_to_file=.true.
	IF (write_to_file) THEN
		IF (i_kpt==1) THEN
			IF (jsp==1) OPEN (679,file='bands_sc_old.1',status='unknown') !This is kind of my birthday 6 july 1992 (S.R.)
			IF (jsp==2) OPEN (680,file='bands_sc_old.2',status='unknown')
		END IF
233
	END IF
234 235 236 237 238

!		write(*,*) 'real zmat size dim 1:', size(zMat%data_r,1), 'dim2:', size(zMat%data_r,2)
!		write(*,*) 'smat dim1', size(smat_unfold%data_r,1), 'dim2', size(smat_unfold%data_r,2),'data',smat_unfold%data_r(2,2)
!		write(222,'(234f15.8)') zMat%data_r
!		write(223,'(234f15.8)') smat_unfold%data_r
239 240


241 242 243 244 245 246 247 248 249 250 251 252 253 254
	IF (zmat%l_real) THEN	
		ALLOCATE(w_n(zMat%matsize2))
	        w_n = 0
!	    IF (method_rubel) THEN
		ALLOCATE(w_n_sum(zMat%matsize2))
	        w_n_sum = 0
!	    END IF
	ELSE
		ALLOCATE(w_n_c(zMat%matsize2))
		w_n_c=0	
!	    IF (method_rubel) THEN
		ALLOCATE(w_n_c_sum(zMat%matsize2))
		w_n_c_sum=0	
!	    END IF
255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271
	END IF
!---------create zmat_s--- smat*zmat---------------------
	select type(zMat)
		type is (t_mat)
		allocate(t_mat::zMat_s)
		select type(zMat_s)
             		type is (t_mat)
	     		zMat_s=zMat
		end select
		type is (t_mpimat)
		allocate(t_mpimat::zMat_s)
		select type(zMat_s)
             		type is (t_mpimat)
	     		zMat_s=zMat
		end select
	end select
!---------------------------------------------------------
272
!		write(345,'(3I6)') lapw%gvec(:,:,jsp)
S.Rost's avatar
S.Rost committed
273
	write (*,*)results%ef
274
        write (*,*) i_kpt
275 276 277 278 279 280 281 282 283 284
	IF (.not. method_rubel) THEN
!          IF (mpi%n_size==1) THEN       
!             call smat_unfold%multiply(zMat,zMat_s)
!          ELSE
!             call smat_unfold%mpimat_multiply(zMat,zMat_s)
!          ENDIF  
           call smat_unfold%multiply(zMat,zMat_s)
        END IF
       !$omp parallel private(j,n_i,nn,na,lo,nk,nki,gi)
       !$omp do
285 286
	DO i=1,zMat%matsize2
		IF (method_rubel) THEN
287
			DO j=1,lapw%nv(jsp)
288 289
				IF (zmat%l_real) THEN
					w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
290
!						write(*,*) 'zMat is real'
291 292
				ELSE
					w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
293
!						write(*,*) 'zMat is complex'
294 295 296 297
				END IF
				IF ((modulo(lapw%gvec(1,j,jsp)+NINT(kpts%sc_list(7,i_kpt)),banddos%s_cell_x)==0).AND.&
				     &(modulo(lapw%gvec(2,j,jsp)+NINT(kpts%sc_list(8,i_kpt)),banddos%s_cell_y)==0).AND.&
				     &(modulo(lapw%gvec(3,j,jsp)+NINT(kpts%sc_list(9,i_kpt)),banddos%s_cell_z)==0)) THEN
298 299
					IF (zmat%l_real) THEN
						w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
300
!							write(*,*) 'zMat is real'
301 302
					ELSE
						w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
303
!							write(*,*) 'zMat is complex'
304
					END IF
305
			   	END IF
306
			END DO
S.Rost's avatar
S.Rost committed
307
!------------------LO's------------------------
308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333
			na=0
			DO n_i=1,atoms%ntype
				DO nn=1,atoms%neq(n_i)
					na=na+1
					DO lo=1,atoms%nlo(n_i)
						nk=lapw%nkvec(lo,na)
						DO nki=1,nk
							gi=lapw%kvec(nki,lo,na)
							j=lapw%nv(jsp)+lapw%index_lo(lo,na)+nki
							IF (zmat%l_real) THEN
								w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
							ELSE
								w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
							END IF
							IF ((modulo(lapw%gvec(1,gi,jsp)+NINT(kpts%sc_list(7,i_kpt)),banddos%s_cell_x)==0).AND.&
							   &(modulo(lapw%gvec(2,gi,jsp)+NINT(kpts%sc_list(8,i_kpt)),banddos%s_cell_y)==0).AND.&
							   &(modulo(lapw%gvec(3,gi,jsp)+NINT(kpts%sc_list(9,i_kpt)),banddos%s_cell_z)==0)) THEN
								IF (zmat%l_real) THEN
									w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
								ELSE
									w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
								END IF
							END IF
						END DO
					END DO
				END DO
S.Rost's avatar
S.Rost committed
334 335
			END DO
!--------------------------LO's finished----------------
336
		ELSE
337
			DO j=1,lapw%nv(jsp)
338
!				DO k=1,zMat%matsize1
339
					IF (zmat%l_real) THEN
340 341
!						w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
						w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat_s%data_r(j,i)
342
					ELSE
343 344
!						w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(k,i)*smat_unfold%data_c(j,k)
						w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat_s%data_c(j,i)
345
					END IF
346
!				END DO
347 348 349
				IF ((modulo(lapw%gvec(1,j,jsp)+NINT(kpts%sc_list(7,i_kpt)),banddos%s_cell_x)==0).AND.&
				   &(modulo(lapw%gvec(2,j,jsp)+NINT(kpts%sc_list(8,i_kpt)),banddos%s_cell_y)==0).AND.&
				   &(modulo(lapw%gvec(3,j,jsp)+NINT(kpts%sc_list(9,i_kpt)),banddos%s_cell_z)==0)) THEN
350
!					DO k=1,zMat%matsize1
351
						IF (zmat%l_real) THEN
352 353
!							w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
							w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat_s%data_r(j,i)
354
						ELSE
355 356
!							w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(k,i)*smat_unfold%data_c(j,k)
							w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat_s%data_c(j,i)
357
						END IF
358
!					END DO
359
				END IF
S.Rost's avatar
S.Rost committed
360
			END DO
361
!			write(1250+mpi%irank,'(4f15.8)') w_n_c(i),w_n_c_sum(i)
S.Rost's avatar
S.Rost committed
362
!------------------LO's------------------------
363 364 365 366 367 368 369 370 371
      			na=0
      			DO n_i=1,atoms%ntype
        			DO nn=1,atoms%neq(n_i)
          				na=na+1
          				DO lo=1,atoms%nlo(n_i)
						nk=lapw%nkvec(lo,na)
						DO nki=1,nk
							gi=lapw%kvec(nki,lo,na)
							j=lapw%nv(jsp)+lapw%index_lo(lo,na)+nki
372
				!			DO k=1,zMat%matsize1
373
								IF (zmat%l_real) THEN
374 375
				!						w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
									w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat_s%data_r(j,i)
376
								ELSE
377 378
				!						w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(k,i)*smat_unfold%data_c(j,k)
									w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat_s%data_c(j,i)
379
								END IF
380
				!			END DO
381 382 383
							IF ((modulo(lapw%gvec(1,gi,jsp)+NINT(kpts%sc_list(7,i_kpt)),banddos%s_cell_x)==0).AND.&
							   &(modulo(lapw%gvec(2,gi,jsp)+NINT(kpts%sc_list(8,i_kpt)),banddos%s_cell_y)==0).AND.&
							   &(modulo(lapw%gvec(3,gi,jsp)+NINT(kpts%sc_list(9,i_kpt)),banddos%s_cell_z)==0)) THEN
384
			!					DO k=1,zMat%matsize1
385
									IF (zmat%l_real) THEN
386 387
			!							w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
										w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat_s%data_r(j,i)
388
									ELSE
389 390
			!							w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(k,i)*smat_unfold%data_c(j,k)
										w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat_s%data_c(j,i)
391
									END IF
392
		    	!					END DO
393 394 395 396
							END IF
						END DO
					END DO
				END DO
397
			END DO
S.Rost's avatar
S.Rost committed
398
!--------------------------LO's finished----------------
399 400
		END IF
!		IF (method_rubel) THEN
401
		IF (write_to_file) THEN
402
			IF (zmat%l_real) THEN
403
				IF (w_n(i)/w_n_sum(i)<0) w_n(i)=0   ! delete negative entries
404 405
				IF (jsp==1) write(679,'(3f15.8)') kpts%sc_list(10,i_kpt), ((eig(i)-results%ef)*hartree_to_ev_const),w_n(i)/w_n_sum(i)
				IF (jsp==2) write(680,'(3f15.8)') kpts%sc_list(10,i_kpt), ((eig(i)-results%ef)*hartree_to_ev_const),w_n(i)/w_n_sum(i)
406
				IF ((w_n(i)/w_n_sum(i)>1).or.(w_n(i)/w_n_sum(i)<0)) write(*,*) 'w_n/sum larger 1 or smaller 0', w_n(i)/w_n_sum(i), 'eigenvalue',eig(i)
407
			ELSE
408
				IF (real(w_n_c(i))<0) w_n_c(i)=0    ! delete negative entries
409 410
				IF (jsp==1) write(679,'(4f15.8)') kpts%sc_list(10,i_kpt), ((eig(i)-results%ef)*hartree_to_ev_const),w_n_c(i)/w_n_c_sum(i)
				IF (jsp==2) write(680,'(4f15.8)') kpts%sc_list(10,i_kpt), ((eig(i)-results%ef)*hartree_to_ev_const),w_n_c(i)/w_n_c_sum(i)
411 412
				IF ((abs(w_n_c(i)/w_n_c_sum(i))>1).or.(real(w_n_c(i))<0)) write(*,*) 'w_n_c/sum larger 1 or smaller 0', w_n_c(i)/w_n_c_sum(i), 'eigenvalue',eig(i)
		        END IF
413 414 415 416 417 418 419 420 421 422
                END IF
		IF (zmat%l_real) THEN
			IF (w_n(i)/w_n_sum(i)<0) w_n(i)=0   ! delete negative entries
			unfoldingBuffer(i,i_kpt,jsp)=w_n(i)/w_n_sum(i)
			IF ((w_n(i)/w_n_sum(i)>1).or.(w_n(i)/w_n_sum(i)<0)) write(*,*) 'w_n/sum larger 1 or smaller 0', w_n(i)/w_n_sum(i), 'eigenvalue',eig(i)
		ELSE
			IF (real(w_n_c(i))<0) w_n_c(i)=0    ! delete negative entries
			unfoldingBuffer(i,i_kpt,jsp)=w_n_c(i)/w_n_c_sum(i)
			IF ((abs(w_n_c(i)/w_n_c_sum(i))>1).or.(real(w_n_c(i))<0)) write(*,*) 'w_n_c/sum larger 1 or smaller 0', w_n_c(i)/w_n_c_sum(i), 'eigenvalue',eig(i)
	        END IF
423 424 425 426 427 428 429 430 431 432 433 434
!		ELSE
!			IF (zmat%l_real) THEN
!				IF (jsp==1) write(679,'(3f15.8)') kpt_dist, ((eig(i)-results%ef)*hartree_to_ev_const),w_n(i)
!				IF (jsp==2) write(680,'(3f15.8)') kpt_dist, ((eig(i)-results%ef)*hartree_to_ev_const),w_n(i)
!				IF ((w_n(i)>1).or.(w_n(i)<0)) write(*,*) 'w_n larger 1 or smaller 0', w_n(i), 'eigenvalue',eig(i)
!			ELSE
!				IF (jsp==1) write(679,'(4f15.8)') kpt_dist, ((eig(i)-results%ef)*hartree_to_ev_const),w_n_c(i)
!				IF (jsp==2) write(680,'(4f15.8)') kpt_dist, ((eig(i)-results%ef)*hartree_to_ev_const),w_n_c(i)
!				IF ((abs(w_n_c(i))>1).or.(real(w_n_c(i))<0)) write(*,*) 'w_n_c larger 1 or smaller 0', w_n_c(i), 'eigenvalue',eig(i)
!	        	END IF
!		END IF			
	END DO
435 436 437
       !$omp end do
       !$omp end parallel
        write (*,*) 'finished',i_kpt
438
	IF (i_kpt==kpts%nkpt) THEN
439
		IF (write_to_file .AND. jsp==1) CLOSE (679)
440
		IF (jsp==input%jspins) THEN
441
			IF (write_to_file .AND. jsp==2) CLOSE (680)
442 443 444 445
			!kpts%bk(:,:)=kpts%sc_list(11:13,:)
			write(*,*) 'Unfolded Bandstructure calculated succesfully, calledby=calculate_plot_w_n'
			!CALL juDFT_error('Unfolded Bandstructure created succesfully - use band_sc.gnu to plot', calledby='calculate_plot_w_n')

446
		END IF
447 448
	END IF
 END SUBROUTINE
449 450 451 452 453 454 455 456 457 458 459 460 461

	
SUBROUTINE write_band_sc(kpts,results,eFermiPrev)
     USE m_types
     USE m_juDFT
     USE m_constants
     IMPLICIT NONE
	TYPE(t_results),INTENT(IN)  :: results
	TYPE(t_kpts),INTENT(IN)     :: kpts
        REAL, INTENT(IN) :: eFermiPrev
	INTEGER :: i,i_kpt,jsp

	OPEN (679,file='bands_sc.1',status='unknown') !This is kind of my birthday 6 july 1992 (S.R.)
462
	IF (SIZE(results%unfolding_weights,3)==2) OPEN (680,file='bands_sc.2',status='unknown')
463 464 465 466 467 468 469 470 471
        DO jsp=1,SIZE(results%unfolding_weights,3)
		DO i_kpt=1,SIZE(results%unfolding_weights,2)
			DO i=1,results%neig(i_kpt,jsp)
				IF (jsp==1) write(679,'(4f15.8)') kpts%sc_list(10,i_kpt), ((results%eig(i,i_kpt,1)-eFermiPrev)*hartree_to_ev_const),results%unfolding_weights(i,i_kpt,1)
				IF (jsp==2) write(680,'(4f15.8)') kpts%sc_list(10,i_kpt), ((results%eig(i,i_kpt,2)-eFermiPrev)*hartree_to_ev_const),results%unfolding_weights(i,i_kpt,2)
			END DO
		END DO
	END DO
	CLOSE (679)
472
	IF (SIZE(results%unfolding_weights,3)==2) CLOSE (680)
473 474
	write(*,*) 'Unfolded Bandstructure written succesfully - use band_sc.gnu to plot, calledby=write_band_sc',eFermiPrev
END SUBROUTINE
475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547
      	
      SUBROUTINE write_gnu_sc(nosyp,d,ssy,input)
      	USE m_types
	USE m_juDFT
      IMPLICIT NONE

      TYPE(t_input),INTENT(IN) :: input
      INTEGER, INTENT (IN) :: nosyp
      REAL,    INTENT (IN) :: d(nosyp)
      CHARACTER(len=1), INTENT (IN) :: ssy(nosyp)
      
      INTEGER n,aoff,adel
      CHARACTER(LEN=200) tempTitle
      aoff = iachar('a')-1
      adel = iachar('a')-iachar('A')
      !write(*,*) aoff,adel 

      OPEN (27,file='band_sc.gnu',status='unknown')
      WRITE (27,*) 'reset'
      WRITE (27,900)
      WRITE (27,901)
      WRITE (27,902)
      WRITE (27,903)
      WRITE(tempTitle,'(10a)') input%comment
      IF(TRIM(ADJUSTL(tempTitle)).EQ.'') THEN
         tempTitle = "Fleur Bandstructure"
      END IF
      WRITE (27,904) TRIM(ADJUSTL(tempTitle))
      DO n = 1, nosyp
        WRITE (27,905) d(n),d(n)
      ENDDO
      WRITE (27,906) d(1),d(nosyp)
!
! nomal labels
!
      IF (iachar(ssy(1)) < aoff ) THEN
        WRITE (27,907) ssy(1),d(1),achar(92)
      ELSE
        WRITE (27,907) " ",d(1),achar(92)
      ENDIF
      DO n = 2, nosyp-1
        IF (iachar(ssy(n)) < aoff ) THEN 
          WRITE (27,908) ssy(n),d(n),achar(92)
        ELSE
          WRITE (27,908) " ",d(n),achar(92)
        ENDIF
      ENDDO
      IF (iachar(ssy(nosyp)) < aoff ) THEN
        WRITE (27,909) ssy(nosyp),d(nosyp)
      ELSE
        WRITE (27,909) " ",d(nosyp)
      ENDIF
!
! greek labels
!
      DO n = 1, nosyp
        IF (iachar(ssy(n)) > aoff ) THEN
          WRITE (27,914) achar(iachar(ssy(n))-adel),d(n)
        ENDIF
      ENDDO
!
! now write the rest
!
      WRITE (27,910)
      WRITE (27,*) 'set palette model RGB'
      WRITE (27,*) 'set palette defined (-2 "black", -1 "white" ,0 "white",',achar(92)
      WRITE (27,*) '0.67 "light-blue",1 "blue")'
      WRITE (27,*) 'set cbrange [-2:1]'
      WRITE (27,*) 'unset colorbox'
      WRITE (27,*) 'size1(x)=0.9*x**(0.4)'
      WRITE (27,*) 'color1(x)=0.3+x/2.4'
      WRITE (27,*) 'size2(x)=0.35*(1-x**(0.01))'
      WRITE (27,*) 'color2(x)=1.15*(x-1)'
548
      WRITE (27,*) 'e_f=0.000000 #fermi energy is already corrected when using hdf5'
549
      WRITE (27,911) d(nosyp)+0.00001,achar(92)
550 551 552 553
      IF (input%jspins == 2) THEN
	WRITE (27,912) achar(92)
	WRITE (27,916) achar(92)
      END IF
554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569
      WRITE (27,913) achar(92)
      WRITE (27,915)
      CLOSE (27)

 900  FORMAT ('set terminal postscript enhanced color "Times-Roman" 20')
 901  FORMAT ('set xlabel ""')
 902  FORMAT ('set ylabel "E - E_F (eV)"')
 903  FORMAT ('set nokey')
 904  FORMAT ('set title "',a,'"')
 905  FORMAT ('set arrow from',f9.5,', -9.0 to',f9.5,',  5.0 nohead')
 906  FORMAT ('set arrow from',f9.5,', 0.0 to',f9.5,', 0.0 nohead lt 3')
 907  FORMAT ('set xtics ("',a1,'"',f9.5,', ',a)
 908  FORMAT ('           "',a1,'"',f9.5,', ',a)
 909  FORMAT ('           "',a1,'"',f9.5,'  )')
 910  FORMAT ('set ytics -8,2,4')
 911  FORMAT ('plot [0:',f9.5,'] [-9:5] ',a)
570 571 572 573
 912  FORMAT ('"bands_sc.2" using 1:($2-e_f):(size1($3)):(color1($3))  w p pt 7 ps variable lc palette, ',a)
 916  FORMAT ('"bands_sc.2" using 1:($2-e_f):(size2($3)):(color2($3)) w p pt 7 ps variable lc palette,',a)
 913  FORMAT ('"bands_sc.1" using 1:($2-e_f):(size1($3)):(color1($3))  w p pt 7 ps variable lc palette, ',a)
 915  FORMAT ('"bands_sc.1" using 1:($2-e_f):(size2($3)):(color2($3)) w p pt 7 ps variable lc palette')
574 575
 914  FORMAT ('set label "',a1,'" at ',f9.5,', -9.65 center font "Symbol,20"')
      END SUBROUTINE write_gnu_sc
576
END MODULE m_unfold_band_kpts