TABLE OF CONTENTS

ABINIT/gwls_LanczosResolvents [ Modules ]

[ Top ] [ Modules ]

NAME

 gwls_LanczosResolvents

FUNCTION

  .

COPYRIGHT

 Copyright (C) 2009-2018 ABINIT group (JLJ, BR, MC)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 
28 module gwls_LanczosResolvents
29 !----------------------------------------------------------------------------------------------------
30 ! This module contains routines to compute the matrix elements of the resolent, namely
31 !
32 !                       < phi_l | 1/ (H-z) | psi_l' >
33 !
34 ! where H is the Hamiltonian, and z should be such that we are far from its poles.
35 !
36 !----------------------------------------------------------------------------------------------------
37 ! local modules
38 
39 use m_gwls_utility
40 use gwls_wf
41 use gwls_TimingLog
42 use gwls_hamiltonian
43 use gwls_lineqsolver
44 use gwls_GWlanczos
45 
46 use defs_basis
47 use m_profiling_abi
48 use m_xmpi
49 use m_io_tools,  only : get_unit
50 
51 implicit none
52 save
53 private

m_hamiltonian/build_preconditioned_Hamiltonian_Lanczos_basis [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  build_preconditioned_Hamiltonian_Lanczos_basis

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_LanczosResolvents

CHILDREN

      zgetrf,zgetri

SOURCE

292 subroutine build_preconditioned_Hamiltonian_Lanczos_basis(seed_vector)
293 !----------------------------------------------------------------------------------------------------
294 ! This function Computes the Lanczos basis of the preconditioned Hamiltonian.
295 !----------------------------------------------------------------------------------------------------
296 
297 !This section has been created automatically by the script Abilint (TD).
298 !Do not modify the following lines by hand.
299 #undef ABI_FUNC
300 #define ABI_FUNC 'build_preconditioned_Hamiltonian_Lanczos_basis'
301 !End of the abilint section
302 
303 implicit none
304 complex(dpc), intent(in) :: seed_vector(npw_g)
305 
306 integer :: l
307 integer :: ierr
308 integer :: mpi_communicator
309 
310 ! *************************************************************************
311 
312 
313 
314 mpi_communicator = mpi_enreg%comm_fft
315 
316 ! Compute the Lanczos basis as well as the eigenvalues of the T matrix
317 ! Seed must also be preconditioned!
318 LR_seeds(:,1) = precondition_one_on_C(:)*seed_vector(:)
319 
320 call block_lanczos_algorithm(mpi_communicator, matrix_function_preconditioned_Hamiltonian, LR_kmax, LR_nseeds, npw_g, LR_seeds, &
321 &                            LR_alpha, LR_beta, Hamiltonian_Qk)
322 
323 
324 call diagonalize_lanczos_banded(LR_kmax,LR_nseeds,npw_g, LR_alpha,LR_beta,  & 
325 Hamiltonian_Qk, LR_Hamiltonian_eigenvalues,.false.)
326 
327 
328 
329 ! Update the Lanczos vectors with the preconditioning matrix
330 do l = 1, LR_kmax
331 Hamiltonian_Qk(:,l) =  precondition_C(:)*Hamiltonian_Qk(:,l)
332 end do
333 
334 ! compute the M matrix
335 
336 ! Compute Q^dagger . Q
337 call ZGEMM('C','N',LR_kmax,LR_kmax,npw_g,cmplx_1,Hamiltonian_Qk,npw_g,Hamiltonian_Qk,npw_g,cmplx_0,LR_M_matrix,LR_kmax)
338 call xmpi_sum(LR_M_matrix,mpi_communicator,ierr) ! sum on all processors working on FFT!
339 
340 do l = 1, LR_kmax
341 LR_M_matrix(l,:) = LR_M_matrix(l,:)*LR_Hamiltonian_eigenvalues(l)
342 end do
343 
344 end subroutine build_preconditioned_Hamiltonian_Lanczos_basis

m_hamiltonian/cleanup_LanczosResolvents [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  cleanup_LanczosResolvents

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_Projected_AT,gwls_Projected_BT

CHILDREN

      zgetrf,zgetri

SOURCE

180 subroutine cleanup_LanczosResolvents
181 !----------------------------------------------------------------------------------------------------
182 !
183 ! This subroutine cleans up global arrays.
184 !
185 !
186 !----------------------------------------------------------------------------------------------------
187 
188 !This section has been created automatically by the script Abilint (TD).
189 !Do not modify the following lines by hand.
190 #undef ABI_FUNC
191 #define ABI_FUNC 'cleanup_LanczosResolvents'
192 !End of the abilint section
193 
194 implicit none
195 
196 ! *************************************************************************
197 
198 ABI_DEALLOCATE(precondition_C)
199 ABI_DEALLOCATE(precondition_one_on_C)
200 ABI_DEALLOCATE(LR_alpha)
201 ABI_DEALLOCATE(LR_beta )
202 ABI_DEALLOCATE(LR_seeds)
203 ABI_DEALLOCATE(Hamiltonian_Qk)
204 ABI_DEALLOCATE(LR_Hamiltonian_eigenvalues)
205 ABI_DEALLOCATE(LR_M_matrix)
206 
207 end subroutine cleanup_LanczosResolvents

m_hamiltonian/compute_resolvent_column_shift_lanczos [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  compute_resolvent_column_shift_lanczos

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_Projected_AT

CHILDREN

      zgetrf,zgetri

SOURCE

367 subroutine compute_resolvent_column_shift_lanczos(nz, list_z, nvec, list_left_vectors, seed_vector, matrix_elements_resolvent)
368 !----------------------------------------------------------------------------------------------------
369 ! This function Computes one column of the resolvent using the shift lanczos algorithm,
370 ! namely
371 !
372 !
373 !               I_l(z) = < left_vector_l | xi(z) > 
374 !
375 ! where |xi(z) > is obtained from the shift lanczos scheme. 
376 !----------------------------------------------------------------------------------------------------
377 
378 !This section has been created automatically by the script Abilint (TD).
379 !Do not modify the following lines by hand.
380 #undef ABI_FUNC
381 #define ABI_FUNC 'compute_resolvent_column_shift_lanczos'
382 !End of the abilint section
383 
384 implicit none
385 
386 integer     , intent(in) :: nz
387 complex(dpc), intent(in) :: list_z(nz)
388 
389 integer     , intent(in) :: nvec
390 complex(dpc), intent(in) :: list_left_vectors(npw_g,nvec)
391 
392 complex(dpc), intent(in) :: seed_vector(npw_g)
393 
394 complex(dpc), intent(out) :: matrix_elements_resolvent(nz,nvec)
395 
396 ! local variables
397 integer        :: iz, l
398 complex(dpc)   :: z
399 integer        :: ierr
400 integer        :: mpi_communicator
401 
402 complex(dpc)   :: right_vec(LR_kmax)
403 complex(dpc)   :: left_vecs(LR_kmax,nvec)
404 complex(dpc)   :: work_vec(LR_kmax)
405 complex(dpc)   :: shift_lanczos_matrix(LR_kmax,LR_kmax)
406 
407 complex(dpc)   :: work_array(npw_g)
408 
409 ! *************************************************************************
410 
411 
412 ! processes will communicate along FFT rows
413 mpi_communicator = mpi_enreg%comm_fft
414 
415 !----------------------------------------------------------------------------------------------------
416 ! First, build the Lanczos basis for the preconditioned Hamiltonian
417 !----------------------------------------------------------------------------------------------------
418 call build_preconditioned_Hamiltonian_Lanczos_basis(seed_vector)
419 
420 
421 !----------------------------------------------------------------------------------------------------
422 ! Next, generate arrays which do not depend on z, the external shift.
423 !----------------------------------------------------------------------------------------------------
424 
425 ! Q^dagger . C^{-2} | seed > 
426 work_array(:) = precondition_one_on_C(:)**2*seed_vector 
427 
428 call ZGEMV('C',npw_g,LR_kmax,cmplx_1,Hamiltonian_Qk,npw_g,work_array,1,cmplx_0,right_vec,1)
429 call xmpi_sum(right_vec,mpi_communicator,ierr) ! sum on all processors working on FFT!
430 
431 !  Q^dagger | left_vectors > 
432 call ZGEMM('C','N',LR_kmax,nvec,npw_g,cmplx_1,Hamiltonian_Qk,npw_g,list_left_vectors,npw_g,cmplx_0,left_vecs,LR_kmax)
433 call xmpi_sum(left_vecs,mpi_communicator,ierr) ! sum on all processors working on FFT!
434 
435 !----------------------------------------------------------------------------------------------------
436 ! Use shift Lanczos to compute all matrix elements!
437 !----------------------------------------------------------------------------------------------------
438 
439 do iz = 1, nz
440 
441 z = list_z(iz)
442 
443 ! Generate the matrix to be inverted
444 shift_lanczos_matrix(:,:) = LR_M_matrix(:,:)
445 
446 do l = 1, LR_kmax
447 shift_lanczos_matrix(l,l) = shift_lanczos_matrix(l,l)-z
448 end do 
449 
450 ! since z could be complex, the matrix is not necessarily hermitian. Invert using general Lapack scheme
451 call invert_general_matrix(LR_kmax,shift_lanczos_matrix)
452 ! the matrix now contains the inverse!
453 
454 
455 ! | work_vec > = M^{-1} . | right_vec >
456 call ZGEMV('N',LR_kmax,LR_kmax,cmplx_1,shift_lanczos_matrix,LR_kmax,right_vec,1,cmplx_0,work_vec,1)
457 
458 
459 ! matrix_elements = < right_vecs | work_vec > 
460 call ZGEMV('C',LR_kmax,nvec,cmplx_1,left_vecs,LR_kmax,work_vec,1,cmplx_0,matrix_elements_resolvent(iz,:),1)
461 
462 end do
463 
464 
465 end subroutine compute_resolvent_column_shift_lanczos

m_hamiltonian/compute_resolvent_column_shift_lanczos_right_vectors [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  compute_resolvent_column_shift_lanczos_right_vectors

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_Projected_BT

CHILDREN

      zgetrf,zgetri

SOURCE

488 subroutine compute_resolvent_column_shift_lanczos_right_vectors(seed_vector, right_vec)
489 !----------------------------------------------------------------------------------------------------
490 ! The goal of this routine is to participate in the computation of a column of the resolvent 
491 ! using the shift lanczos algorithm. This routine should be called FIRST.
492 !
493 ! This routine:
494 !
495 !       1) builds the preconditioned Hamiltonian Lanczos basis, using the seed_vector, and 
496 !          stores result in the Module variables.
497 !
498 !       2) prepares and returns the Lanczos basis-projected right vectors, ready for further processing.
499 !
500 !
501 !----------------------------------------------------------------------------------------------------
502 
503 !This section has been created automatically by the script Abilint (TD).
504 !Do not modify the following lines by hand.
505 #undef ABI_FUNC
506 #define ABI_FUNC 'compute_resolvent_column_shift_lanczos_right_vectors'
507 !End of the abilint section
508 
509 implicit none
510 
511 complex(dpc), intent(in) :: seed_vector(npw_g)
512 complex(dpc), intent(out):: right_vec(LR_kmax)
513 
514 ! local variables
515 integer        :: mpi_communicator
516 integer        :: ierr
517 
518 complex(dpc)   :: work_array(npw_g)
519 
520 ! *************************************************************************
521 
522 
523 ! processes will communicate along FFT rows
524 mpi_communicator = mpi_enreg%comm_fft
525 
526 !----------------------------------------------------------------------------------------------------
527 ! First, build the Lanczos basis for the preconditioned Hamiltonian
528 !----------------------------------------------------------------------------------------------------
529 call build_preconditioned_Hamiltonian_Lanczos_basis(seed_vector)
530 
531 !----------------------------------------------------------------------------------------------------
532 ! Next, generate arrays which do not depend on z, the external shift.
533 !----------------------------------------------------------------------------------------------------
534 
535 ! Q^dagger . C^{-2} | seed > 
536 work_array(:) = precondition_one_on_C(:)**2*seed_vector(:)
537 
538 call ZGEMV('C', npw_g, LR_kmax, cmplx_1, Hamiltonian_Qk, npw_g, work_array, 1, cmplx_0, right_vec, 1)
539 call xmpi_sum(right_vec,mpi_communicator,ierr) ! sum on all processors working on FFT!
540 
541 
542 end subroutine compute_resolvent_column_shift_lanczos_right_vectors

m_hamiltonian/invert_general_matrix [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  invert_general_matrix

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_LanczosResolvents,gwls_Projected_BT

CHILDREN

      zgetrf,zgetri

SOURCE

566 subroutine invert_general_matrix(n,matrix)
567 !----------------------------------------------------------------------------------------------------
568 ! Simple wrapper around lapack routines to invert a general matrix
569 !----------------------------------------------------------------------------------------------------
570 
571 !This section has been created automatically by the script Abilint (TD).
572 !Do not modify the following lines by hand.
573 #undef ABI_FUNC
574 #define ABI_FUNC 'invert_general_matrix'
575 !End of the abilint section
576 
577 implicit none
578 
579 integer,      intent(in)    :: n
580 complex(dpc), intent(inout) :: matrix(n,n)
581 
582 integer :: info
583 integer :: ipiv(n)
584 
585 integer        :: debug_unit
586 character(50)  :: debug_filename
587 
588 
589 
590 
591 complex(dpc) :: work(n)
592 
593 ! *************************************************************************
594 
595 call ZGETRF( n, n, matrix, n, ipiv, info )
596 if ( info /= 0) then        
597   debug_unit = get_unit()
598   write(debug_filename,'(A,I4.4,A)') 'LAPACK_DEBUG_PROC=',mpi_enreg%me,'.log'
599 
600   open(debug_unit,file=trim(debug_filename),status='unknown')
601 
602   write(debug_unit,'(A)')      '*********************************************************************************************'
603   write(debug_unit,'(A,I4,A)') '*      ERROR: info = ',info,' in ZGETRF(1), gwls_LanczosResolvents'
604   write(debug_unit,'(A)')      '*********************************************************************************************'
605 
606   close(debug_unit)
607 
608 end if
609 
610 
611 call ZGETRI( n, matrix, n, ipiv, work, n, info )
612 
613 if ( info /= 0) then        
614   debug_unit = get_unit()
615   write(debug_filename,'(A,I4.4,A)') 'LAPACK_DEBUG_PROC=',mpi_enreg%me,'.log'
616 
617   open(debug_unit,file=trim(debug_filename),status='unknown')
618 
619   write(debug_unit,'(A)')      '*********************************************************************************************'
620   write(debug_unit,'(A,I4,A)') '*      ERROR: info = ',info,' in ZGETRI(1), gwls_LanczosResolvents'
621   write(debug_unit,'(A)')      '*********************************************************************************************'
622 
623   close(debug_unit)
624 
625 end if
626 
627 
628 
629 
630 end subroutine invert_general_matrix

m_hamiltonian/matrix_function_preconditioned_Hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  matrix_function_preconditioned_Hamiltonian

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

CHILDREN

      zgetrf,zgetri

SOURCE

228 subroutine matrix_function_preconditioned_Hamiltonian(vector_out,vector_in,Hsize)
229 !----------------------------------------------------------------------------------------------------
230 ! This subroutine is a simple wrapper around the preconditioned Hamiltonian operator which acts as
231 !
232 !                               C^{-1} . H . C^{-1}
233 !
234 ! where C^{-2} ~ 1 / T, where T is the kinetic energy.
235 !
236 !----------------------------------------------------------------------------------------------------
237 
238 !This section has been created automatically by the script Abilint (TD).
239 !Do not modify the following lines by hand.
240 #undef ABI_FUNC
241 #define ABI_FUNC 'matrix_function_preconditioned_Hamiltonian'
242 !End of the abilint section
243 
244 implicit none
245 integer,      intent(in)  :: Hsize
246 complex(dpc), intent(out) :: vector_out(Hsize)
247 complex(dpc), intent(in)  :: vector_in(Hsize)
248 
249 ! local variables
250 real(dp)  :: psikg(2,npw_g)
251 
252 complex(dpc)  :: tmp_vector(Hsize)
253 
254 ! *************************************************************************
255 
256 ! convert from one format to the other
257 
258 tmp_vector(:) = precondition_one_on_C(:)*vector_in(:)
259 
260 psikg(1,:) = dble (tmp_vector(:))
261 psikg(2,:) = dimag(tmp_vector(:))
262 
263 call Hpsik(psikg)
264 
265 tmp_vector(:) = cmplx_1*psikg(1,:)+cmplx_i*psikg(2,:)
266 
267 vector_out(:) = precondition_one_on_C(:)*tmp_vector(:)
268 
269 
270 end subroutine matrix_function_preconditioned_Hamiltonian

m_hamiltonian/setup_LanczosResolvents [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  setup_LanczosResolvents

FUNCTION

  .

INPUTS

OUTPUT

PARENTS

      gwls_Projected_AT,gwls_Projected_BT

CHILDREN

      zgetrf,zgetri

SOURCE

103 subroutine setup_LanczosResolvents(kmax, prec)
104 !----------------------------------------------------------------------------------------------------
105 !
106 ! This subroutine prepares global work arrays in order to use the Lanczos scheme to 
107 ! compute matrix elements of inverse operators.
108 !
109 !
110 !----------------------------------------------------------------------------------------------------
111 
112 !This section has been created automatically by the script Abilint (TD).
113 !Do not modify the following lines by hand.
114 #undef ABI_FUNC
115 #define ABI_FUNC 'setup_LanczosResolvents'
116 !End of the abilint section
117 
118 implicit none
119 
120 !------------------------------
121 ! input/output variables
122 !------------------------------
123 
124 integer, intent(in) :: kmax  ! number of Lanczos steps
125 logical, intent(in) :: prec  ! should there be preconditioning?
126 
127 ! *************************************************************************
128 
129 
130 LR_kmax   = kmax
131 LR_nseeds = 1      ! only one seed
132 
133 ! Prepare the algorithm
134 
135 ABI_ALLOCATE(precondition_C,        (npw_g))
136 ABI_ALLOCATE(precondition_one_on_C, (npw_g))
137 
138 
139 if ( prec ) then
140   call set_precondition()
141 
142   precondition_one_on_C(:) = sqrt(pcon(:))
143   !precondition_one_on_C(:) = pcon(:) ! yields very bad results!
144   precondition_C(:)        = cmplx_1/precondition_one_on_C(:)
145 
146 else 
147   precondition_C(:)        = cmplx_1
148   precondition_one_on_C(:) = cmplx_1
149 end if
150 
151 ABI_ALLOCATE(LR_alpha,       (LR_nseeds,LR_nseeds,LR_kmax))
152 ABI_ALLOCATE(LR_beta ,       (LR_nseeds,LR_nseeds,LR_kmax))
153 ABI_ALLOCATE(LR_seeds,       (npw_g,LR_nseeds))
154 ABI_ALLOCATE(Hamiltonian_Qk, (npw_g,LR_kmax))
155 ABI_ALLOCATE(LR_Hamiltonian_eigenvalues, (LR_kmax))
156 ABI_ALLOCATE(LR_M_matrix, (LR_kmax,LR_kmax))
157 
158 end subroutine setup_LanczosResolvents