Branch data 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 : : !!***
|
