Line data Source code
1 : !!****f* ABINIT/m_rot_cg
2 : !! NAME
3 : !! m_rot_cg
4 : !!
5 : !! FUNCTION
6 : !! Rotate the cg coefficient with the rotation matrix obtained from the
7 : !! diagonalization of the non-diagonal occupation matrix produced by DMFT.
8 : !!
9 : !! INPUTS
10 : !!
11 : !! OUTPUT
12 : !!
13 : !! SOURCE
14 : !!
15 : !! TODO /!\ No parallel computing yet !
16 :
17 : #if defined HAVE_CONFIG_H
18 : #include "config.h"
19 : #endif
20 :
21 : #include "abi_common.h"
22 :
23 : module m_rot_cg
24 :
25 : use defs_basis
26 : use m_abi_linalg, only : abi_xgemm
27 : use m_abicore
28 : use m_errors
29 : use m_xmpi
30 :
31 : implicit none
32 :
33 : private
34 :
35 : public :: rot_cg
36 :
37 : contains
38 : !!***
39 :
40 : !!****f* ABINIT/diag_occ
41 : !! NAME
42 : !! diag_occ
43 : !!
44 : !! FUNCTION
45 : !! Use for DMFT in KGB parallelisation. Diagonalise the occupation matrix
46 : !! and return diagonalised occupations and associated eigenvectors sorted
47 : !! with descending occupation
48 : !!
49 : !! INPUTS
50 : !! occ_nd_cpx(nband,nband) = matrix of non diagonal occupations for DMFT
51 : !! nband = number of bands to be processed
52 : !!
53 : !! OUTPUT
54 : !! occ_diag(2,nband) = diagonal occupations in the new band space
55 : !!
56 : !! SOURCE
57 : !!
58 : !! TODO /!\ No parallel computing yet !
59 : !! TODO add the possibility of using ScaLAPACK to do computation in parallel
60 :
61 220 : subroutine diag_occ(occ_nd_cpx,nband,occ_diag)
62 :
63 : !Arguments ------------------------------------
64 : !scalars
65 : integer, intent(in) :: nband
66 : !! type(MPI_type),intent(inout) :: mpi_enreg
67 : !! type(dataset_type),intent(in) :: dtset
68 : !! type(paw_dmft_type), intent(in) :: paw_dmft
69 : !no_abirules
70 : complex(dp), intent(inout) :: occ_nd_cpx(nband,nband)
71 : real(dp), intent(inout) :: occ_diag(nband)
72 : !Local variables-------------------------------
73 : integer :: info,lwork
74 : character(len=500) :: message
75 440 : real(dp) :: rwork(3*nband-1)
76 220 : complex(dp), allocatable :: work(:)
77 : ! *************************************************************************
78 :
79 : DBG_ENTER("COLL")
80 :
81 : !! Use the opposite to have zheev orders the eigenvalues by descending order.
82 : !! Afterwards, we multiply by -1 once again.
83 34540 : occ_nd_cpx(:,:) = - occ_nd_cpx(:,:)
84 :
85 : !! Get diagonal occupations and associated base
86 :
87 : ! Compute the optimal working array size
88 220 : ABI_MALLOC(work,(1))
89 220 : call zheev('v','u',nband,occ_nd_cpx(:,:),nband,occ_diag(:),work(:),-1,rwork(:),info)
90 220 : lwork = int(work(1))
91 220 : ABI_FREE(work)
92 :
93 : ! Compute the eigenvalues (occ_diag) and vectors
94 660 : ABI_MALLOC(work,(lwork))
95 :
96 220 : call zheev('v','u',nband,occ_nd_cpx(:,:),nband,occ_diag(:),work(:),lwork,rwork(:),info)
97 :
98 : !! Obtain the true eigenvalues of occupation matrix in descending order
99 2860 : occ_diag(:) = - occ_diag(:)
100 :
101 220 : ABI_FREE(work)
102 :
103 220 : if (info > 0) then
104 0 : message = ""
105 : write(message,"(a,i5)") " something wrong happened with the diagonalization of &
106 0 : & the occupation matrix (didn't converge), info=",info
107 0 : ABI_ERROR(message)
108 220 : else if (info < 0) then
109 0 : message = ""
110 : write(message,"(a,i5)") " something wrong happened with the diagonalization of &
111 0 : & the occupation matrix (bad input argument), info=",info
112 0 : ABI_ERROR(message)
113 : end if
114 :
115 : DBG_EXIT("COLL")
116 :
117 220 : end subroutine diag_occ
118 : !!***
119 :
120 : !!****f* ABINIT/rot_cg
121 : !! NAME
122 : !! rot_cg
123 : !!
124 : !! FUNCTION
125 : !! Use for DMFT in KGB parallelisation. Diagonalise the occupation matrix
126 : !! and use the resulting base to represent the wave functions.
127 : !!
128 : !! INPUTS
129 : !! occ_nd(2,nband,nband) = matrix of non diagonal occupations for DMFT
130 : !! cwavef(2,npw,nband) = Fourier coefficients of wave functions for all bands
131 : !! npw = number of G vectors computed in this iteration
132 : !! nband = number of bands to be processed
133 : !! blocksize = size of the block for the LO.. algorithm
134 : !! still has to be equal to nband
135 : !! nspinor = number of spinor components
136 : !! first_bandc = index of the first correlated band
137 : !! nbandc = number of correlated bands
138 : !!
139 : !! OUTPUT
140 : !! occ_diag(nband) = diagonal occupations in the new band space
141 : !!
142 : !! SIDE EFFECT
143 : !! cwavef is rotated with the unitary matrix obtained from the diagonalization
144 : !! of occupations (occ_nd)
145 : !! SOURCE
146 : !!
147 : !! TODO /!\ No parallel computing yet !
148 : !! TODO add the possibility of using ScaLAPACK to do computation in parallel
149 : !! TODO Make the computation of the new wf parallel
150 :
151 220 : subroutine rot_cg(occ_nd,cwavef,npw,nband,blocksize,nspinor,first_bandc,nbandc,occ_diag,dmft_optim)
152 :
153 : !Arguments ------------------------------------
154 : !scalars
155 : integer, intent(in) :: blocksize,first_bandc,nband,nbandc,npw,nspinor
156 : logical, intent(in) :: dmft_optim
157 : !! type(MPI_type),intent(inout) :: mpi_enreg
158 : !! type(dataset_type),intent(in) :: dtset
159 : !! type(paw_dmft_type), intent(in) :: band_in
160 : !no_abirules
161 : real(dp), intent(in) :: occ_nd(2,blocksize,blocksize)
162 : real(dp), intent(inout) :: occ_diag(blocksize)
163 : real(dp), intent(inout) :: cwavef(2,npw,blocksize,nspinor)
164 : !Local variables-------------------------------
165 : !scalars
166 : integer :: ispinor,n
167 : character(len=500) :: message
168 : !arrays
169 220 : real(dp), allocatable :: occ_diag_red(:)
170 220 : complex(dp), allocatable :: mat_tmp(:,:),mat_tmp2(:,:),occ_nd_cpx(:,:)
171 : !complex(kind=dp) :: cwavef_rot_g(nbandc, nspinor)
172 : ! *************************************************************************
173 :
174 : DBG_ENTER("COLL")
175 :
176 220 : if (nband /= blocksize) then
177 0 : message = " DMFT in KGB cannot be used with multiple blocks yet. Make sure that bandpp*npband = nband."
178 0 : ABI_ERROR(message)
179 : end if
180 :
181 : !! Initialization
182 :
183 880 : ABI_MALLOC(mat_tmp,(npw,nbandc))
184 660 : ABI_MALLOC(mat_tmp2,(npw,nbandc))
185 660 : ABI_MALLOC(occ_diag_red,(nbandc))
186 880 : ABI_MALLOC(occ_nd_cpx,(nbandc,nbandc))
187 :
188 : occ_nd_cpx(:,:) = cmplx(occ_nd(1,first_bandc:first_bandc+nbandc-1,first_bandc:first_bandc+nbandc-1), &
189 34540 : & occ_nd(2,first_bandc:first_bandc+nbandc-1,first_bandc:first_bandc+nbandc-1),kind=dp)
190 :
191 : !! Get diagonal occupations and associated base
192 :
193 220 : call diag_occ(occ_nd_cpx(:,:),nbandc,occ_diag_red(:))
194 :
195 7260 : do n=1,nband
196 7260 : if (n < first_bandc .or. n >= first_bandc+nbandc) then
197 4400 : occ_diag(n) = occ_nd(1,n,n)
198 : else
199 2640 : occ_diag(n) = occ_diag_red(n-first_bandc+1)
200 : end if
201 : end do ! n
202 :
203 : !! Compute the corresponding wave functions if nothing wrong happened
204 : ! $c^{rot}_{n,k}(g) = \sum_{n'} [\bar{f_{n',n}} * c_{n',k}(g)]$
205 :
206 : ! Correct a bug in the formula when using TRIQS
207 220 : if (dmft_optim) occ_nd_cpx(:,:) = conjg(occ_nd_cpx(:,:))
208 :
209 440 : do ispinor=1,nspinor
210 : mat_tmp(:,:) = cmplx(cwavef(1,1:npw,first_bandc:first_bandc+nbandc-1,ispinor), &
211 497332 : & cwavef(2,1:npw,first_bandc:first_bandc+nbandc-1,ispinor),kind=dp)
212 : call abi_xgemm("n","n",npw,nbandc,nbandc,cone,mat_tmp(:,:),npw, &
213 220 : & occ_nd_cpx(:,:),nbandc,czero,mat_tmp2(:,:),npw)
214 497332 : cwavef(1,1:npw,first_bandc:first_bandc+nbandc-1,ispinor) = dble(mat_tmp2(:,:))
215 497552 : cwavef(2,1:npw,first_bandc:first_bandc+nbandc-1,ispinor) = aimag(mat_tmp2(:,:))
216 : end do ! ispinor
217 :
218 : !do ig=1,npw
219 : ! cwavef_rot_g(:,:) = czero
220 : ! do n=1,nbandc
221 : ! do np=1,nbandc
222 : ! cwavef_rot_g(n,:) = cwavef_rot_g(n,:) + occ_nd_cpx(np, n) * &
223 : !& cmplx(cwavef(1,ig,np+first_bandc-1,:), cwavef(2,ig,np+first_bandc-1,:), kind=dp)
224 : ! end do
225 : ! end do
226 : ! cwavef(1,ig,first_bandc:first_bandc+nbandc-1,:) = dreal(cwavef_rot_g)
227 : ! cwavef(2,ig,first_bandc:first_bandc+nbandc-1,:) = dimag(cwavef_rot_g)
228 : ! end do
229 :
230 220 : ABI_FREE(mat_tmp)
231 220 : ABI_FREE(mat_tmp2)
232 220 : ABI_FREE(occ_diag_red)
233 220 : ABI_FREE(occ_nd_cpx)
234 :
235 : DBG_EXIT("COLL")
236 :
237 220 : end subroutine rot_cg
238 : !!***
239 :
240 : end module m_rot_cg
241 : !!***
|
