1! Copyright (C) 2006-2008 Dmitry Korotin dmitry@korotin.name
2! This file is distributed under the terms of the
3! GNU General Public License. See the file `License'
4! in the root directory of the present distribution,
5! or http://www.gnu.org/copyleft/gpl.txt .
6!
7#define ZERO (0.d0,0.d0)
8#define ONE (1.d0,0.d0)
9
10!-----------------------------------------------------------------------
11SUBROUTINE write_hamiltonian_default(nwan,hamk,iunhamilt)
12!-----------------------------------------------------------------------
13
14  USE io_global, ONLY: stdout
15  USE kinds, ONLY: DP
16  USE constants,  ONLY : rytoev
17  USE klist, ONLY: nks, wk
18
19  IMPLICIT NONE
20  INTEGER, INTENT(in) :: nwan, iunhamilt
21  COMPLEX(DP), INTENT(in) :: hamk(nwan,nwan,nks)
22
23  INTEGER :: i,j, ik, ios
24  COMPLEX(DP), ALLOCATABLE :: hamk2(:,:)
25  real(DP) :: eps = 1.d-8, hr,hi
26
27  WRITE(stdout,'(/5x,a36,i5,a9)') 'Hamiltonian is in the default format,', nks, 'k-points'
28  WRITE(stdout,'(5x,a48/)') 'ATTENTION: All k-points weights are real weights'
29
30  OPEN (iunhamilt, file = 'hamilt', status = 'unknown', form = 'formatted', err = 300, iostat = ios)
31300 CALL errore ('HMLT', 'Opening hamilt', abs (ios) )
32
33  WRITE(iunhamilt,*) nks,nwan
34
35  ALLOCATE(hamk2(nwan,nwan))
36
37  DO ik = 1, nks
38
39     WRITE(iunhamilt,'(f15.12)') wk(ik)
40
41     hamk2 = hamk(:,:,ik) * rytoev
42
43     DO i=1, nwan
44        DO j=1, nwan
45
46           hr = abs(dreal(hamk2(i,j)))
47           hi = abs(dimag(hamk2(i,j)))
48           IF((hr>=eps).and.(hi>=eps)) WRITE(iunhamilt,'(2f13.8)') dreal(hamk2(i,j)), aimag(hamk2(i,j))
49           IF ((hr<eps).and.(hi>=eps)) WRITE(iunhamilt,'(f3.0,f13.8)') 0., aimag(hamk2(i,j))
50           IF ((hr>=eps).and.(hi<eps)) WRITE(iunhamilt,'(f13.8,f3.0)') dreal(hamk2(i,j)), 0.
51           IF ((hr<eps).and.(hi<eps)) WRITE(iunhamilt,'(2f3.0)') 0., 0.
52
53        ENDDO
54     ENDDO
55
56  ENDDO
57
58  DEALLOCATE(hamk2)
59
60  CLOSE(iunhamilt)
61
62END SUBROUTINE write_hamiltonian_default
63
64!-----------------------------------------------------------------------
65SUBROUTINE write_hamiltonian_amulet(nwan,hamk,hash,iunhamilt)
66!-----------------------------------------------------------------------
67! Special hamiltonian format for the AMULET code instegration
68! www.amulet-code.org
69
70  USE io_global, ONLY: stdout
71  USE kinds, ONLY: DP
72  USE constants,  ONLY : rytoev
73  USE klist, ONLY: nks, wk, xk
74  USE lsda_mod, ONLY : nspin
75  USE run_info, ONLY : title
76  USE global_version, ONLY : version_number
77
78  IMPLICIT NONE
79  INTEGER, INTENT(in) :: nwan, hash, iunhamilt
80  COMPLEX(DP), INTENT(in) :: hamk(nwan,nwan,nks)
81
82  INTEGER :: i,j, ik, ios
83  COMPLEX(DP), ALLOCATABLE :: hamk2(:,:)
84  real(DP) :: eps = 1.d-8, hr,hi
85  CHARACTER(LEN=9)  :: cdate, ctime
86
87  CALL date_and_tim( cdate, ctime )
88
89  WRITE(stdout,'(/5x,a36,i5,a9)') 'Hamiltonian is in the AMULET format,', nks/nspin, 'k-points'
90  WRITE(stdout,'(5x,a48/)') 'ATTENTION: All k-points weights are real weights'
91
92  OPEN (iunhamilt, file = 'hamilt.am', status = 'unknown', form = 'formatted', err = 300, iostat = ios)
93300 CALL errore ('HMLT', 'Opening hamilt', abs (ios) )
94
95  write(iunhamilt,'(a30,2a10/)') '# This file was generated on: ', cdate, ctime
96  IF( trim(title) .NE. '' ) write(iunhamilt,'(a2,a80/)') '# ', title
97
98  WRITE(iunhamilt,'(a10)') '&codestamp'
99  WRITE(iunhamilt,'(a3,a6)') 'QE_', version_number
100  WRITE(iunhamilt,*)
101
102  WRITE(iunhamilt,'(a5)') '&hash'
103  WRITE(iunhamilt,*) hash
104  WRITE(iunhamilt,*)
105
106  WRITE(iunhamilt,'(a6)') '&nspin'
107  WRITE(iunhamilt,'(i1)') nspin
108  WRITE(iunhamilt,*)
109
110  WRITE(iunhamilt,'(a4)') '&nkp'
111  WRITE(iunhamilt,'(i5)') nks/nspin
112  WRITE(iunhamilt,*)
113
114  WRITE(iunhamilt,'(a4)') '&dim'
115  WRITE(iunhamilt,'(i3)') nwan
116  WRITE(iunhamilt,*)
117
118  WRITE(iunhamilt,'(a8)') '&kpoints'
119  DO ik=1, nks/nspin
120    WRITE(iunhamilt,'(f15.12,3f9.5)') wk(ik), xk(:,ik)
121  END DO
122  WRITE(iunhamilt,*)
123
124  WRITE(iunhamilt,'(a12)') '&hamiltonian'
125
126  ALLOCATE(hamk2(nwan,nwan))
127
128  DO ik = 1, nks
129
130     hamk2 = hamk(:,:,ik) * rytoev
131
132     DO i=1, nwan
133        DO j=i, nwan
134
135           hr = abs(dreal(hamk2(i,j)))
136           hi = abs(dimag(hamk2(i,j)))
137           IF((hr>=eps).and.(hi>=eps)) WRITE(iunhamilt,'(2f13.8)') dreal(hamk2(i,j)), aimag(hamk2(i,j))
138           IF ((hr<eps).and.(hi>=eps)) WRITE(iunhamilt,'(f3.0,f13.8)') 0., aimag(hamk2(i,j))
139           IF ((hr>=eps).and.(hi<eps)) WRITE(iunhamilt,'(f13.8,f3.0)') dreal(hamk2(i,j)), 0.
140           IF ((hr<eps).and.(hi<eps)) WRITE(iunhamilt,'(2f3.0)') 0., 0.
141
142        ENDDO
143     ENDDO
144
145  ENDDO
146
147  DEALLOCATE(hamk2)
148
149  CLOSE(iunhamilt)
150
151END SUBROUTINE write_hamiltonian_amulet
152
153!-----------------------------------------------------------------------
154SUBROUTINE write_systemdata_amulet(hash,nelec,iunsystem)
155!-----------------------------------------------------------------------
156! Damp of the system data for the AMULET code instegration
157! www.amulet-code.org
158
159  USE io_global, ONLY: stdout
160  USE kinds, ONLY: DP
161  USE constants,  ONLY : rytoev
162  USE klist, ONLY: nks, wk, xk
163  USE lsda_mod, ONLY : nspin
164  USE run_info, ONLY : title
165  USE ions_base, ONLY : nat, atm, tau, ityp
166  USE cell_base, ONLY : alat, at
167  USE ener, ONLY : ef
168  USE wannier_new, ONLY : nwan, wan_in
169  USE global_version, ONLY : version_number
170
171  IMPLICIT NONE
172  INTEGER, INTENT(in) :: hash,iunsystem
173  REAL(DP), INTENT(in) :: nelec
174  INTEGER :: ios, i, j
175  CHARACTER(LEN=9)  :: cdate, ctime
176  CHARACTER :: l_symb(4)
177  DATA l_symb/'s','p','d','f'/
178  INTEGER :: orbitals(7,4)
179  DATA orbitals/ 1, 0, 0, 0, 0, 0, 0, &
180                 3, 4, 2, 0, 0, 0, 0, &
181                 7, 8, 6, 9, 5, 0, 0, &
182                 13, 14, 12, 15, 11, 16, 10 /
183  INTEGER, PARAMETER :: nwannierblocksmax = 25
184
185  INTEGER :: nblocks, block_dim(nwannierblocksmax), block_l(nwannierblocksmax), &
186             block_atom(nwannierblocksmax), block_wannier(nwannierblocksmax,9), block_start(nwannierblocksmax)
187
188  interface
189    SUBROUTINE split_basis_into_blocks(nblocks,block_dim,block_l,block_atom,block_wannier,block_start)
190      INTEGER, INTENT(OUT) :: nblocks, block_dim(:), block_atom(:), block_l(:), block_wannier(:,:), block_start(:)
191    END SUBROUTINE split_basis_into_blocks
192  end interface
193
194  CALL date_and_tim( cdate, ctime )
195
196  OPEN (iunsystem, file = 'system.am', status = 'unknown', form = 'formatted', err = 300, iostat = ios)
197300 CALL errore ('HMLT', 'Opening system.am', abs (ios) )
198
199  write(iunsystem,'(a30,2a10/)') '# This file was generated on: ', cdate, ctime
200  IF( trim(title) .NE. '' ) write(iunsystem,'(a2,a80/)') '# ', title
201
202  WRITE(iunsystem,'(a5)') '&hash'
203  WRITE(iunsystem,*) hash
204  WRITE(iunsystem,*)
205
206  WRITE(iunsystem,'(a10)') '&codestamp'
207  WRITE(iunsystem,'(a3,a6)') 'QE_',version_number
208  WRITE(iunsystem,*)
209
210  WRITE(iunsystem,'(a5)') '&cell'
211  WRITE(iunsystem,'(f12.9)') alat
212  DO i=1,3
213    WRITE(iunsystem,'(3f9.5)') at(:,i)
214  END DO
215  WRITE(iunsystem,*)
216
217  WRITE(iunsystem,'(a6)') '&atoms'
218  WRITE(iunsystem,'(i5)') nat
219  DO i=1, nat
220    WRITE(iunsystem,'(a4,1x,3f9.5)') atm(ityp(i)), tau(:,i)
221  END DO
222  WRITE(iunsystem,*)
223
224  WRITE(iunsystem,'(a6)') '&nelec'
225  WRITE(iunsystem,'(f7.2)') nelec
226  WRITE(iunsystem,*)
227
228  WRITE(iunsystem,'(a7)') '&efermi'
229  WRITE(iunsystem,'(f8.4)') ef*rytoev
230  WRITE(iunsystem,*)
231
232  CALL split_basis_into_blocks(nblocks,block_dim,block_l,block_atom,block_wannier,block_start)
233  WRITE(iunsystem,'(a20)') '# Basis description:'
234  WRITE(iunsystem,'(a14)') '# dim, nblocks'
235  WRITE(iunsystem,'(a74)') '# atom_sym, atom_num, l_sym, block_dim, block_start, orbitals(1:block_dim)'
236  WRITE(iunsystem,'(a6)') '&basis'
237  WRITE(iunsystem,'(i2,i4)') nwan, nblocks
238  DO i=1, nblocks
239    WRITE(iunsystem,'(a3,i3,a2,i2,i4,4x,7i2)') atm(ityp(block_atom(i))), block_atom(i), l_symb(block_l(i)+1), &
240                       block_dim(i), block_start(i), ( orbitals(wan_in(j,1)%ing(1)%m,block_l(i)+1), &
241                        j=block_wannier(i,1), block_wannier(i,block_dim(i)) )
242  END DO
243  WRITE(iunsystem,*)
244
245  CLOSE(iunsystem)
246
247END SUBROUTINE write_systemdata_amulet
248
249SUBROUTINE split_basis_into_blocks(nblocks,block_dim,block_l,block_atom,block_wannier,block_start)
250  USE kinds, only : DP
251  USE wannier_new, only : nwan, wan_in
252  USE io_global, ONLY: stdout
253
254  implicit none
255  INTEGER :: i,j, iblock
256  INTEGER, INTENT(OUT) :: nblocks, block_dim(:), block_atom(:), block_l(:), block_wannier(:,:), block_start(:)
257
258  nblocks = 0
259  block_dim = 0
260  block_l = -1
261  block_atom = 0
262  block_wannier = 0
263
264  iblock = 1
265  block_start(iblock) = 1
266
267  j=1
268  DO i = 1, nwan-1
269    block_wannier(iblock,j) = i
270    j = j+1
271    IF( ((wan_in(i,1)%iatom .ne. wan_in(i+1,1)%iatom) .OR. wan_in(i,1)%ing(1)%l .ne. wan_in(i+1,1)%ing(1)%l) &
272        ) THEN
273      block_dim(iblock) = i - block_start(iblock) + 1
274      block_atom(iblock) = wan_in(i,1)%iatom
275      block_l(iblock) = wan_in(i,1)%ing(1)%l
276      iblock = iblock + 1
277      block_start(iblock) = i+1
278      j=1
279    END IF
280  END DO
281
282  ! the last wannier
283  block_dim(iblock) = nwan - block_start(iblock) + 1
284  block_atom(iblock) = wan_in(nwan,1)%iatom
285  block_l(iblock) = wan_in(nwan,1)%ing(1)%l
286  block_wannier(iblock,j) = i
287
288  nblocks = iblock
289
290!   write(stdout,'(/5x,a19)') 'Blocks of orbitals:'
291!   DO iblock = 1, nblocks
292!     write(stdout,*) 'Block', iblock, block_dim(iblock), block_atom(iblock), &
293!                             block_l(iblock), block_wannier(iblock,1:block_dim(iblock))
294!   END DO
295
296END SUBROUTINE split_basis_into_blocks
297