Line data Source code
1 : !!****m* ABINIT/m_bseinterp
2 : !! NAME
3 : !! m_bseinterp
4 : !!
5 : !! FUNCTION
6 : !!
7 : !! COPYRIGHT
8 : !! Copyright (C) 2014-2026 ABINIT group (Y. Gillet, M.Giantomassi)
9 : !! This file is distributed under the terms of the
10 : !! GNU General Public License, see ~abinit/COPYING
11 : !! or http://www.gnu.org/copyleft/gpl.txt .
12 : !!
13 : !! SOURCE
14 :
15 : #if defined HAVE_CONFIG_H
16 : #include "config.h"
17 : #endif
18 :
19 : #include "abi_common.h"
20 :
21 : MODULE m_bseinterp
22 :
23 : use defs_basis
24 : use m_abicore
25 : use m_bs_defs
26 : use m_xmpi
27 : use m_errors
28 : use m_nctk
29 : use m_haydock_io
30 : use m_linalg_interfaces
31 : use netcdf
32 :
33 : use m_fstrings, only : indent, strcat, sjoin, itoa
34 : use defs_datatypes, only : pseudopotential_type
35 : use m_hide_blas, only : xdotc
36 : use m_fft_mesh, only : calc_ceigr
37 : use m_crystal, only : crystal_t
38 : use m_bz_mesh, only : kmesh_t
39 : use m_double_grid, only : double_grid_t, get_kpt_from_indices_coarse
40 : use m_wfd, only : wfdgw_t
41 : use m_pawtab, only : pawtab_type
42 :
43 : implicit none
44 :
45 : private
46 : !!***
47 :
48 : !----------------------------------------------------------------------
49 :
50 : !!****t* m_haydock/interpolator_t
51 : !! NAME
52 : !! interpolator_t
53 : !!
54 : !! FUNCTION
55 : !! Store the overlap matrix elements needed for the interpolation of the BSE Hamiltonian
56 : !!
57 : !! TODO
58 : !! Decide if we want to make the number of bands k-dependent.
59 : !!
60 : !! SOURCE
61 :
62 : type,public :: interpolator_t
63 :
64 : integer :: nvert=8
65 : ! Number of vertices for interpolation
66 :
67 : integer :: method
68 : ! Interpolation method (YG or Rohlfing & Louie or ...)
69 :
70 : integer :: mband_dense, mband_coarse
71 : ! Max number of bands dense and coarse
72 :
73 : integer :: nsppol
74 : ! Number of spin channels
75 :
76 : integer, allocatable :: corresp(:,:,:)
77 : ! corresp(max_nreh,nvert,spin)
78 : ! it_coarse, idiv -> it_coarse (idiv-th neighbour)
79 :
80 : real(dp),allocatable :: interp_factors(:,:)
81 : ! interp_factors(nvert,ndiv)
82 : ! index_in_fine_box -> k-point in Trans_interp
83 :
84 : complex(gwp),allocatable :: overlaps(:,:,:,:,:)
85 : ! Overlaps between dense and coarse mesh
86 : ! overlaps(mband_coarse,mband_dense,ivertex_coarse,double_grid%nkpt_dense,spin)
87 :
88 : complex(dp),allocatable :: btemp(:), ctemp(:)
89 : ! Temporary arrays for work
90 :
91 : ! Pointers to datatypes that are already in memory
92 : type(double_grid_t),pointer :: double_grid => null()
93 : ! Mapping between coarse and dense mesh
94 :
95 : contains
96 : procedure :: init => interpolator_init ! Construct the object
97 : procedure :: free => interpolator_free ! Free memory
98 : procedure :: normalize => interpolator_normalize ! Normalize the overlaps
99 : procedure :: alloc_work => int_alloc_work ! Alloc temp memory
100 : procedure :: int_free => int_free ! Free temp memory
101 :
102 : end type interpolator_t
103 : !!***
104 :
105 : !----------------------------------------------------------------------
106 :
107 : CONTAINS !=======================================================================
108 : !!***
109 :
110 : !!****f* m_bseinterp/interpolator_init
111 : !! NAME
112 : !! interpolator_init
113 : !!
114 : !! FUNCTION
115 : !! Construct the interpolator object
116 : !!
117 : !! INPUTS
118 : !!
119 : !! OUTPUT
120 : !!
121 : !! SOURCE
122 :
123 4 : subroutine interpolator_init(interpolator, double_grid, Wfd_dense, Wfd_coarse, &
124 4 : Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab, method)
125 :
126 : !Arguments ---------------------------
127 : !scalars
128 : class(interpolator_t),intent(inout) :: interpolator
129 : integer,intent(in) :: method
130 : type(double_grid_t),intent(in),target :: double_grid
131 : type(wfdgw_t),intent(inout) :: Wfd_dense, Wfd_coarse
132 : type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse
133 : type(excparam),intent(in) :: BSp
134 : type(crystal_t),intent(in) :: Cryst
135 : type(pseudopotential_type),intent(in) :: Psps
136 : !arrays
137 : type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw)
138 :
139 : !Local variables ---------------------
140 : !scalars
141 : integer :: nsppol, nvert
142 : integer :: maxnreh, nreh1, nreh2
143 : integer :: mbandc, mbandd, nbzd
144 : real(dp),parameter :: threshold = 0.1_dp
145 : !arrays
146 : character(len=500) :: msg
147 :
148 : !*****************************************************************************
149 :
150 4 : ABI_CHECK(Wfd_coarse%usepaw==0, "PAW not yet supported")
151 4 : ABI_CHECK(BSp%nsppol==1, "nsppol != 1 not yet implemented")
152 4 : ABI_CHECK(Wfd_coarse%nspinor==1, "nspinor != 1 not supported")
153 :
154 : ABI_UNUSED(Pawtab(1)%basis_size)
155 : !paw_overlap(cprj1,cprj2,typat,pawtab,spinor_comm) result(onsite)
156 :
157 4 : interpolator%double_grid => double_grid
158 4 : interpolator%mband_dense = Wfd_dense%mband
159 4 : interpolator%mband_coarse = Wfd_coarse%mband
160 4 : interpolator%method = method
161 4 : interpolator%nsppol = BSp%nsppol
162 :
163 : SELECT CASE(method)
164 : CASE (BSE_INTERP_YG)
165 : nvert = 8
166 : CASE (BSE_INTERP_RL2)
167 : nvert = 2
168 : CASE (BSE_INTERP_RL)
169 0 : nvert = 1
170 : CASE DEFAULT
171 0 : write(msg,'(a,i0)') "Wrong interpolation method: ",method
172 0 : ABI_ERROR(msg)
173 : END SELECT
174 :
175 4 : interpolator%nvert = nvert
176 :
177 4 : mbandc = interpolator%mband_coarse
178 4 : mbandd = interpolator%mband_dense
179 4 : nbzd = double_grid%nbz_dense
180 4 : nsppol = interpolator%nsppol
181 28 : ABI_MALLOC(interpolator%overlaps,(mbandc,mbandd,nvert,nbzd,nsppol))
182 :
183 : call int_compute_overlaps(interpolator,double_grid, Wfd_dense, Wfd_coarse, Kmesh_dense, &
184 4 : & Kmesh_coarse, BSp, Cryst, Psps, Pawtab)
185 :
186 16 : ABI_MALLOC(interpolator%interp_factors,(nvert,double_grid%ndiv))
187 :
188 4 : call int_preprocess_tables(interpolator,double_grid)
189 :
190 4 : nreh1 = BSp%nreh(1)
191 4 : nreh2 = nreh1; if(BSp%nsppol == 2) nreh2 = BSp%nreh(2)
192 4 : maxnreh = MAX(nreh1,nreh2)
193 :
194 20 : ABI_MALLOC(interpolator%corresp,(maxnreh,interpolator%nvert,interpolator%nsppol))
195 :
196 4 : call int_compute_corresp(interpolator,BSp,double_grid)
197 :
198 4 : end subroutine interpolator_init
199 : !!***
200 :
201 : !-------------------------------------------------------------------
202 :
203 : !!****f* m_bseinterp/int_alloc_work
204 : !! NAME
205 : !! int_alloc_work
206 : !!
207 : !! FUNCTION
208 : !! Allocate temporary arrays
209 : !!
210 : !! INPUTS
211 : !!
212 : !! SOURCE
213 :
214 2 : subroutine int_alloc_work(interpolator, work_size)
215 :
216 : !Arguments ---------------------------
217 : !scalars
218 : class(interpolator_t),intent(inout) :: interpolator
219 : integer,intent(in) :: work_size
220 : !*****************************************************************************
221 :
222 6 : ABI_MALLOC(interpolator%btemp,(work_size))
223 4 : ABI_MALLOC(interpolator%ctemp,(work_size))
224 :
225 2 : end subroutine int_alloc_work
226 : !!***
227 :
228 : !-------------------------------------------------------------------
229 :
230 : !!****f* m_bseinterp/int_free
231 : !! NAME
232 : !! int_free
233 : !!
234 : !! FUNCTION
235 : !! Deallocate temporary arrays
236 : !!
237 : !! SOURCE
238 :
239 2 : subroutine int_free(interpolator)
240 :
241 : !Arguments ---------------------------
242 : class(interpolator_t),intent(inout) :: interpolator
243 : !*****************************************************************************
244 :
245 2 : ABI_SFREE(interpolator%btemp)
246 2 : ABI_SFREE(interpolator%ctemp)
247 :
248 2 : end subroutine int_free
249 : !!***
250 :
251 : !-------------------------------------------------------------------
252 :
253 : !!****f* m_bseinterp/int_compute_overlaps
254 : !! NAME
255 : !! int_compute_overlaps
256 : !!
257 : !! FUNCTION
258 : !! Compute the overlaps prefactors
259 : !!
260 : !! INPUTS
261 : !!
262 : !! OUTPUT
263 : !!
264 : !! SOURCE
265 :
266 4 : subroutine int_compute_overlaps(interpolator, double_grid, Wfd_dense, Wfd_coarse, &
267 4 : Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab)
268 :
269 : !Arguments ---------------------------
270 : !scalars
271 : class(interpolator_t),intent(inout) :: interpolator
272 : type(double_grid_t),intent(in),target :: double_grid
273 : type(wfdgw_t),intent(inout) :: Wfd_dense, Wfd_coarse
274 : type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse
275 : type(excparam),intent(in) :: BSp
276 : type(crystal_t),intent(in) :: Cryst
277 : type(pseudopotential_type),intent(in) :: Psps
278 : !arrays
279 : type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw)
280 :
281 : !Local variables ---------------------
282 : !scalars
283 : integer :: nprocs, my_rank, ierr
284 : integer :: nfft, nspinor, nsppol, nvert
285 : integer :: ib_coarse, ib_dense, ik_coarse, ik_dense
286 : integer :: spin, iorder, ivertex, ix, iy, iz, bstart, bstop
287 : real(dp),parameter :: threshold = 0.1_dp
288 : complex(gwp) :: ovlp
289 : !arrays
290 : integer :: curindices_dense(6), curindices_coarse(3)
291 : integer :: neighbour(3)
292 : integer :: g0(3),g01(3),diffg0(3)
293 4 : complex(gwp),allocatable :: ur_coarse(:),ur_dense(:)
294 4 : complex(gwp),allocatable :: ceigr(:)
295 : !*****************************************************************************
296 :
297 4 : nprocs = xmpi_comm_size(Wfd_coarse%comm)
298 4 : my_rank = xmpi_comm_rank(Wfd_coarse%comm)
299 :
300 : ABI_UNUSED(Pawtab(1)%basis_size)
301 :
302 : ! Ensure Wfd and Wfd_coarse use the same FFT mesh.
303 4 : call wfd_dense%change_ngfft(Cryst,Psps,Wfd_coarse%ngfft)
304 4 : nfft = Wfd_coarse%nfft
305 4 : nspinor = Wfd_coarse%nspinor
306 4 : nsppol = Bsp%nsppol
307 4 : nvert = interpolator%nvert
308 :
309 : ! Allocate workspace for wavefunctions in real space.
310 12 : ABI_MALLOC(ur_coarse,(nfft*nspinor))
311 8 : ABI_MALLOC(ur_dense,(nfft*nspinor))
312 8 : ABI_MALLOC(ceigr,(nfft*nspinor))
313 :
314 117064 : interpolator%overlaps = czero
315 :
316 : ! TODO
317 : ! 1) Choose whether we want to compute only dvv, dcc or all dbb
318 : ! 2) Check the ordering of the loops
319 : ! 3) Improve vertex -> neighbour (in double_grid ?)
320 8 : do spin = 1,nsppol
321 264 : do ik_dense = 1,double_grid%nbz_dense
322 :
323 : ! MPI parallelization
324 : ! We assume that each node owns in memory the full set of wavefunctions
325 : ! both coarse and dense k-mesh and both spins.
326 256 : if (mod(ik_dense, nprocs) /= my_rank) cycle
327 :
328 : ! From ik_dense -> indices_dense
329 256 : iorder = double_grid%iktoint_dense(ik_dense)
330 1024 : g01 = double_grid%g0_dense(:,iorder)
331 1792 : curindices_dense = double_grid%indices_dense(:,iorder)
332 :
333 1860 : do ivertex = 1,nvert
334 :
335 : ! From vertex to neighbour
336 : ! TODO improve this part + permit to choose other neighbour (e.g. nearest neighbour for RL)
337 1600 : if(nvert > 1) then
338 1536 : ix = (ivertex-1)/4
339 1536 : iy = (ivertex-ix*4-1)/2
340 1536 : iz = (ivertex-ix*4-iy*2-1)
341 : else
342 64 : ix = (BSp%rl_nb-1)/4
343 64 : iy = (BSp%rl_nb-ix*4-1)/2
344 64 : iz = (BSp%rl_nb-ix*4-iy*2-1)
345 : end if
346 :
347 6400 : neighbour = [ix,iy,iz]
348 :
349 : ! From indices_dense -> indices_coarse
350 6400 : curindices_coarse = curindices_dense(1:3) + neighbour(:)
351 :
352 : ! From indices_coarse -> ik_ibz in the coarse mesh
353 : call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,&
354 1600 : & double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0)
355 :
356 : ! Take into account a possible umklapp between k_dense and k_coarse
357 6400 : diffg0 = g0 - g01
358 :
359 4456 : if (ANY(diffg0/=0)) then
360 : ! WARNING works only with nspinor = 1 !!!
361 888 : call calc_ceigr(diffg0,nfft,nspinor,Wfd_coarse%ngfft,ceigr)
362 : end if
363 :
364 13056 : do ib_dense = BSp%lomo_spin(spin), BSp%humo_spin(spin)
365 : ! ur(ib_dense, ik_dense)
366 11200 : call wfd_dense%sym_ur(Cryst,Kmesh_dense,ib_dense,ik_dense,spin,ur_dense)
367 :
368 31192 : if (ANY(diffg0/=0)) then
369 : !ur_kbz = ur_kbz*e(ig0r)
370 25466952 : ur_dense(:) = ur_dense(:)*ceigr(:)
371 : end if
372 :
373 : ! Uncomment for the complete overlap
374 : !bstart = BSp%lomo_spin(spin); bstop = BSp%humo_spin(spin)
375 :
376 : ! Compute only dvv or dcc
377 11200 : if (ib_dense <= BSp%homo_spin(spin)) then
378 : ! if ib_dense is a valence band => loop on valence bands
379 4800 : bstart = BSp%lomo_spin(spin); bstop = BSp%homo_spin(spin)
380 : else
381 : ! if ib_dense is a conduction band => loop on conduction bands
382 6400 : bstart = BSp%lumo_spin(spin); bstop = BSp%humo_spin(spin)
383 : end if
384 :
385 52800 : do ib_coarse = bstart, bstop
386 : ! ur(ib_coarse, ik_coarse)
387 40000 : call wfd_coarse%sym_ur(Cryst,Kmesh_coarse,ib_coarse,ik_coarse,spin,ur_coarse)
388 :
389 : ! ovlp = < u_{ib_coarse,ik_coarse} | u_{ib_dense,ik_dense} >
390 40000 : ovlp = xdotc(nfft,ur_coarse,1,ur_dense,1)/nfft
391 :
392 : ! Filter too low values
393 40000 : if (ABS(ovlp) < threshold) ovlp = czero
394 :
395 51200 : interpolator%overlaps(ib_coarse,ib_dense,ivertex,ik_dense,spin) = ovlp
396 : end do ! ib_coarse
397 :
398 : !DBYG
399 : ! write(std_out,*) "nb = ",neighbour
400 : ! write(std_out,*) "(i1,i2,i3,j1,j2,j3) = ",curindices_dense
401 : ! write(std_out,*) "ib = ",ib_dense
402 : ! write(std_out,*) "Sum of dbb = ",REAL(SUM(GWPC_CONJG(interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))*interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))); call flush(std_out)
403 : !ENDDBYG
404 :
405 : end do ! ib_dense
406 : end do ! ivertex
407 : end do ! ik_dense
408 : end do ! spin
409 :
410 4 : ABI_FREE(ur_coarse)
411 4 : ABI_FREE(ur_dense)
412 4 : ABI_FREE(ceigr)
413 :
414 : ! Gather results on each node.
415 4 : call xmpi_sum(interpolator%overlaps,Wfd_coarse%comm,ierr)
416 :
417 4 : end subroutine int_compute_overlaps
418 : !!***
419 :
420 : !----------------------------------------------------------------------
421 :
422 : !!****f* m_bseinterp/int_preprocess_tables
423 : !! NAME
424 : !! int_preprocess_tables
425 : !!
426 : !! FUNCTION
427 : !! Pre-process tables to improve interpolation technique
428 : !!
429 : !! INPUTS
430 : !!
431 : !! OUTPUT
432 : !!
433 : !! SOURCE
434 :
435 4 : subroutine int_preprocess_tables(interpolator,double_grid)
436 :
437 : !Argument ------------------------------------
438 : !scalars
439 : class(interpolator_t),intent(inout) :: interpolator
440 : type(double_grid_t),intent(in) :: double_grid
441 :
442 : !Local variables -----------------------------
443 : !scalars
444 : integer :: iorder,ik_dense,ik_coarse
445 : integer :: ix,iy,iz,ineighbour,curdim, curj
446 : real(dp) :: interp_factor
447 : !arrays
448 : integer :: allxyz(3),curindices_dense(6)
449 4 : integer,allocatable :: curindex(:)
450 : !*********************************************
451 :
452 12 : ABI_MALLOC(curindex,(double_grid%nbz_coarse))
453 36 : curindex = 1
454 :
455 236 : interpolator%interp_factors = zero
456 :
457 260 : do ik_dense = 1,double_grid%nbz_dense
458 :
459 : ! From ik_ibz in the dense mesh -> indices_dense
460 256 : iorder = double_grid%iktoint_dense(ik_dense)
461 : !g01 = double_grid%g0_dense(:,iorder)
462 :
463 : ! From indices_dense -> indices_coarse
464 1792 : curindices_dense = double_grid%indices_dense(:,iorder)
465 :
466 256 : ik_coarse = double_grid%dense_to_coarse(ik_dense)
467 :
468 : ! Compute multi-linear interpolation factors
469 : ! Loop over the neighbours
470 1856 : do ineighbour = 1,interpolator%nvert
471 : !TODO helper function from [ix,iy,iz] -> ineighbour and vice versa
472 1600 : ix = (ineighbour-1)/4
473 1600 : iy = (ineighbour-ix*4-1)/2
474 1600 : iz = (ineighbour-ix*4-iy*2-1)
475 6400 : allxyz = [ix,iy,iz]
476 1600 : interp_factor = one
477 6400 : do curdim = 1,3
478 4800 : if (interpolator%method == BSE_INTERP_RL) then
479 : cycle
480 4608 : else if(interpolator%method == BSE_INTERP_RL2) then
481 0 : if (curdim /= 3) cycle
482 : end if
483 4608 : curj = curindices_dense(3+curdim)
484 : interp_factor = interp_factor*((allxyz(curdim)*(curj*1.0/double_grid%kmult(curdim)))&
485 6400 : & +((1-allxyz(curdim))*(1-(curj*1.0/double_grid%kmult(curdim)))))
486 : end do
487 1856 : interpolator%interp_factors(ineighbour,curindex(ik_coarse)) = interp_factor
488 : end do
489 :
490 260 : curindex(ik_coarse) = curindex(ik_coarse) + 1
491 : end do
492 :
493 4 : ABI_FREE(curindex)
494 :
495 4 : end subroutine int_preprocess_tables
496 : !!***
497 :
498 : !-------------------------------------------------------------------
499 :
500 : !!****f* m_haydock/int_compute_corresp
501 : !! NAME
502 : !! int_compute_corresp
503 : !!
504 : !! FUNCTION
505 : !!
506 : !! INPUTS
507 : !! BSp<type(excparam)=The parameter for the Bethe-Salpeter run.
508 : !! grid <double_grid_t>=Correspondence between coarse and fine k-grid
509 : !! spin=Spin index.
510 : !!
511 : !! OUTPUT
512 : !! corresp(Bsp%nreh(spin),8)= Correspondence between a transition on the
513 : !! coarse mesh and its i-th neighbour for i in [1,2,..,8].
514 : !!
515 : !! TODO:
516 : !! Some operations are faster if we allocate with shape (8,nreh(spin))
517 : !!
518 : !! SOURCE
519 :
520 4 : subroutine int_compute_corresp(interpolator,BSp,double_grid)
521 :
522 : !Arguments ------------------------------------
523 : class(interpolator_t),intent(inout) :: interpolator
524 : type(excparam),intent(in) :: BSp
525 : type(double_grid_t),intent(in) :: double_grid
526 :
527 : !Local variables ------------------------------
528 : !scalars
529 : integer :: spin
530 : integer :: itt,ik_dense,ik_coarse,iorder,it_coarse
531 : integer :: ic,iv,ik_coarse0,it_coarse0,iovlp,ix,iy,iz
532 : !arrays
533 : integer :: curindices_dense(6),curindices_coarse(3),g0(3),g01(3),neighbour(3)
534 : !************************************************************************
535 :
536 8 : do spin=1,interpolator%nsppol
537 3080 : do itt=1,BSp%nreh_interp(spin)
538 : ! From dense itt -> ik_dense, ic, iv
539 3072 : ik_dense = BSp%Trans_interp(itt,spin)%k
540 3072 : ic = BSp%Trans_interp(itt,spin)%c
541 3072 : iv = BSp%Trans_interp(itt,spin)%v
542 :
543 : ! From ik_dense -> indices_dense
544 3072 : iorder = double_grid%iktoint_dense(ik_dense)
545 3072 : g01 = double_grid%g0_dense(:,iorder)
546 :
547 : ! Index of the k-point in the coarse mesh.
548 3072 : ik_coarse0 = double_grid%dense_to_coarse(ik_dense)
549 3072 : it_coarse0 = BSp%vcks2t(iv,ic,ik_coarse0,spin)
550 :
551 : ! From indices_dense -> indices_coarse
552 21504 : curindices_dense = double_grid%indices_dense(:,iorder)
553 :
554 : ! Loop over the 8 neighbors.
555 22276 : do iovlp = 1,interpolator%nvert
556 :
557 : !TODO : helper function from [ix,iy,iz] -> iovlp and vice versa
558 19200 : if(interpolator%nvert > 1) then
559 18432 : ix = (iovlp-1)/4
560 18432 : iy = (iovlp-ix*4-1)/2
561 18432 : iz = (iovlp-ix*4-iy*2-1)
562 : else
563 768 : ix = (BSp%rl_nb-1)/4
564 768 : iy = (BSp%rl_nb-ix*4-1)/2
565 768 : iz = (BSp%rl_nb-ix*4-iy*2-1)
566 : end if
567 76800 : neighbour = [ix,iy,iz]
568 :
569 76800 : curindices_coarse = curindices_dense(1:3) + neighbour(:)
570 :
571 : ! From indices_coarse -> ik_ibz in the coarse mesh
572 : call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,&
573 19200 : & double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0)
574 :
575 : ! From ik_coarse, ic, iv to it_coarse
576 19200 : it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin)
577 :
578 22272 : interpolator%corresp(it_coarse0,iovlp,spin) = it_coarse
579 : end do
580 : end do ! itt
581 : end do
582 :
583 4 : end subroutine int_compute_corresp
584 : !!***
585 :
586 : !----------------------------------------------------------------------
587 :
588 : !!****f* m_bseinterp/interpolator_normalize
589 : !! NAME
590 : !! interpolator_normalize
591 : !!
592 : !! FUNCTION
593 : !! Normalize the overlaps so that \sum_{ib} | d_{kk'}^{b,ib} | ^2 = 1
594 : !!
595 : !! INPUTS
596 : !!
597 : !! OUTPUT
598 : !!
599 : !! SOURCE
600 :
601 4 : subroutine interpolator_normalize(interpolator)
602 :
603 : !Arguments ---------------------------
604 : class(interpolator_t),intent(inout) :: interpolator
605 :
606 : !Local variables ---------------------
607 : !scalars
608 : integer :: spin, ivertex, ib_dense, ik_dense
609 : complex(gwp) :: sum_ovlp
610 : !arrays
611 4 : complex(gwp),allocatable :: overlaps(:)
612 : !*****************************************************************************
613 :
614 12 : ABI_MALLOC(overlaps,(interpolator%mband_coarse))
615 8 : do spin = 1, interpolator%nsppol
616 33 : do ivertex = 1, interpolator%nvert
617 229 : do ib_dense = 1, interpolator%mband_dense
618 13025 : do ik_dense = 1, interpolator%double_grid%nbz_dense
619 115200 : overlaps(:) = interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin)
620 115200 : sum_ovlp = SQRT(REAL(SUM(GWPC_CONJG(overlaps(:))*overlaps(:))))
621 12800 : if (ABS(sum_ovlp) > tol6) then
622 99936 : overlaps(:) = overlaps(:)/sum_ovlp
623 : else
624 15264 : overlaps(:) = czero
625 : end if
626 115400 : interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin) = overlaps(:)
627 : end do
628 : end do
629 : end do
630 : end do
631 4 : ABI_FREE(overlaps)
632 :
633 4 : end subroutine interpolator_normalize
634 : !!***
635 :
636 : !-------------------------------------------------------------------
637 :
638 : !!****f* m_bseinterp/interpolator_free
639 : !! NAME
640 : !! interpolator_free
641 : !!
642 : !! FUNCTION
643 : !! Destroy the interpolator object in memory
644 : !!
645 : !! SOURCE
646 :
647 23 : subroutine interpolator_free(interpolator)
648 :
649 : !Arguments ---------------------------
650 : class(interpolator_t),intent(inout) :: interpolator
651 : !*****************************************************************************
652 :
653 23 : ABI_SFREE(interpolator%overlaps)
654 23 : ABI_SFREE(interpolator%corresp)
655 23 : ABI_SFREE(interpolator%interp_factors)
656 23 : if( associated(interpolator%double_grid) ) nullify(interpolator%double_grid)
657 :
658 23 : end subroutine interpolator_free
659 : !!***
660 :
661 : !-------------------------------------------------------------------
662 :
663 0 : end module m_bseinterp
664 : !!***
|
