unfoldBandKPTS.f90 19.2 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 50 51

    p_kpts=kpts
    write(1088,*) 'banddos%unfoldband: ', banddos%unfoldband
    write(1088,*) 'brav. matrix: '
52 53 54
    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)
55
    write(1088,*) 'brav. rez. matrix: '
56 57 58
    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)
59
    write(1088,*) ' primitive brav. matrix: '
60 61 62
    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)
63
    write(1088,*) 'primitive brav. rez. matrix: '
64
    write(89,'(3f15.8)') p_cell%bmat
65 66 67
    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)
68
    write(1088,'(a,i7,a,i7)') 'kpts%nkpt',kpts%nkpt,'   p_kpts%nkpt',p_kpts%nkpt
69
    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(10,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 101 102

    CALL inv3(cell%bmat,rez_inv_to_internal,rez_inv_det)
    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 163 164 165 166 167 168
	!  	kpts%bk(:,i)=matmul(rez_inv_to_internal,pc_kpoint_c)
	    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
169
    END DO
170 171
    write(91,'(3f15.8)') kpts%bk
    write(92,*) kpts%wtkpt
172
    ALLOCATE (kpts%sc_list(10,p_kpts%nkpt))
173
    kpts%sc_list=list
174
    write(90,'(10f15.8)') kpts%sc_list
175 176
  END SUBROUTINE find_supercell_kpts

177
 SUBROUTINE calculate_plot_w_n(banddos,cell,kpts,smat_unfold,zMat,lapw,i_kpt,jsp,eig,results,input,atoms)
178 179 180 181 182 183 184
	USE m_types
	USE m_juDFT
	USE m_inv3
	USE m_types_mpimat
        USE m_constants
	implicit none

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

	CALL build_primitive_cell(banddos,p_cell,cell)

        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
216 217 218 219
	IF (i_kpt==1) THEN
		IF (jsp==1) OPEN (679,file='bands_sc.1',status='unknown') !This is kind of my birthday 6 july 1992 (S.R.)
		IF (jsp==2) OPEN (680,file='bands_sc.2',status='unknown')
	END IF
220 221 222 223 224

!		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
225

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

228 229 230 231 232 233 234 235 236 237 238 239 240 241 242
	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
	END IF	
243
!		write(345,'(3I6)') lapw%gvec(:,:,jsp)
S.Rost's avatar
S.Rost committed
244
	write (*,*)results%ef
245 246 247
        write (*,*) i_kpt
	DO i=1,zMat%matsize2
		IF (method_rubel) THEN
248
			DO j=1,lapw%nv(jsp)
249 250
				IF (zmat%l_real) THEN
					w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
251
!						write(*,*) 'zMat is real'
252 253
				ELSE
					w_n_c_sum(i)=w_n_c_sum(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
254
!						write(*,*) 'zMat is complex'
255 256 257 258
				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
259 260
					IF (zmat%l_real) THEN
						w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(j,i)
261
!							write(*,*) 'zMat is real'
262 263
					ELSE
						w_n_c(i)=w_n_c(i)+CONJG(zMat%data_c(j,i))*zMat%data_c(j,i)
264
!							write(*,*) 'zMat is complex'
265
					END IF
266
			   	END IF
267
			END DO
S.Rost's avatar
S.Rost committed
268
!------------------LO's------------------------
269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294
			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
295 296
			END DO
!--------------------------LO's finished----------------
297
		ELSE
298
			DO j=1,lapw%nv(jsp)
S.Rost's avatar
S.Rost committed
299
				DO k=1,zMat%matsize1
300
					IF (zmat%l_real) THEN
301
						w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
302
					ELSE
303
						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)
304 305
					END IF
				END DO
306 307 308 309 310 311 312 313 314 315 316
				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
					DO k=1,zMat%matsize1
						IF (zmat%l_real) THEN
							w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
						ELSE
							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)
						END IF
					END DO
				END IF
S.Rost's avatar
S.Rost committed
317 318
			END DO
!------------------LO's------------------------
319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348
      			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
							DO k=1,zMat%matsize1
								IF (zmat%l_real) THEN
									w_n_sum(i)=w_n_sum(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
								ELSE
									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)
								END IF
							END DO
							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
								DO k=1,zMat%matsize1
									IF (zmat%l_real) THEN
										w_n(i)=w_n(i)+zMat%data_r(j,i)*zMat%data_r(k,i)*smat_unfold%data_r(j,k)
									ELSE
										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)
									END IF
		    						END DO
							END IF
						END DO
					END DO
				END DO
349
			END DO
S.Rost's avatar
S.Rost committed
350
!--------------------------LO's finished----------------
351 352
		END IF
!		IF (method_rubel) THEN
353
			IF (zmat%l_real) THEN
354
				IF (w_n(i)/w_n_sum(i)<0) w_n(i)=0   ! delete negative entries
355 356
				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)
357
				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)
358
			ELSE
359
				IF (real(w_n_c(i))<0) w_n_c(i)=0    ! delete negative entries
360 361
				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)
362 363
				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
364 365 366 367 368 369 370 371 372 373 374 375
!		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
376
	IF (i_kpt==kpts%nkpt) THEN
377 378 379 380 381
		IF (jsp==1) CLOSE (679)
		IF (jsp==input%jspins) THEN
			IF (jsp==2) CLOSE (680)
			CALL juDFT_error('Unfolded Bandstructure created succesfully - use band_sc.gnu to plot', calledby='calculate_plot_w_n')
		END IF
382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456
	END IF
 END SUBROUTINE
      	
      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)'
457
      WRITE (27,*) 'e_f=0.000000 #fermi energy is already corrected when using hdf5'
458
      WRITE (27,911) d(nosyp)+0.00001,achar(92)
459 460 461 462
      IF (input%jspins == 2) THEN
	WRITE (27,912) achar(92)
	WRITE (27,916) achar(92)
      END IF
463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478
      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)
479 480 481 482
 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')
483 484
 914  FORMAT ('set label "',a1,'" at ',f9.5,', -9.65 center font "Symbol,20"')
      END SUBROUTINE write_gnu_sc
485
END MODULE m_unfold_band_kpts