TABLE OF CONTENTS


ABINIT/m_pawfgrtab [ Modules ]

[ Top ] [ Modules ]

NAME

  m_pawfgrtab

FUNCTION

  This module contains the definition of the pawfgrtab_type structured datatype,
  as well as related functions and methods.
  pawfgrtab_type variables contain various data expressed on the fine grid
  for a given atom.

COPYRIGHT

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

NOTES

  FOR DEVELOPPERS: in order to preserve the portability of libPAW library,
  please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt

SOURCE

24 #include "libpaw.h"
25 
26 MODULE m_pawfgrtab
27 
28  USE_DEFS
29  USE_MSG_HANDLING
30  USE_MPI_WRAPPERS
31  USE_MEMORY_PROFILING
32 
33  use m_paral_atom, only : get_my_atmtab, free_my_atmtab, get_my_natom
34 
35  implicit none
36 
37  private

m_pawfgrtab/pawfgrtab_copy [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

  pawfgrtab_copy

FUNCTION

  Copy one pawfgrtab datastructure into another
  Can take into accound changes of dimensions
  Can copy a shared pawfgrtab into distributed ones (when parallelism is activated)

INPUTS

  mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc
  comm_atom=--optional-- MPI communicator over atoms
  pawfgrtab_in(:)<type(pawfgrtab_type)>= input paw_an datastructure

SIDE EFFECTS

  pawfgrtab_copy(:)<type(pawfgrtab_type)>= output pawfgrtab datastructure

NOTES

  paw_fgrtab_copy must have been allocated in the calling function.

PARENTS

      m_pawfgrtab

CHILDREN

SOURCE

435 subroutine pawfgrtab_copy(pawfgrtab_in,pawfgrtab_cp, &
436 &                         mpi_atmtab,comm_atom) ! optional arguments
437 
438 
439 !This section has been created automatically by the script Abilint (TD).
440 !Do not modify the following lines by hand.
441 #undef ABI_FUNC
442 #define ABI_FUNC 'pawfgrtab_copy'
443 !End of the abilint section
444 
445  implicit none
446 
447 !Arguments ------------------------------------
448 !scalars
449 integer,optional,intent(in) :: comm_atom
450 !arrays
451  integer,optional,target,intent(in) :: mpi_atmtab(:)
452  type(Pawfgrtab_type),intent(in) :: pawfgrtab_in(:)
453  type(Pawfgrtab_type),intent(inout),target :: pawfgrtab_cp(:)
454 
455 !Local variables-------------------------------
456 !scalars
457  integer :: ij,ij1,istat,l_size,l_size2,my_comm_atom,my_natom,nfgd,nspden,paral_case
458  integer :: siz_in, siz_max, siz_out
459  logical ::  my_atmtab_allocated,paral_atom
460  character(len=500) :: msg
461 !arrays
462  integer,pointer :: my_atmtab(:)
463  type(Pawfgrtab_type), pointer :: pawfgrtab_out(:)
464 
465 ! *************************************************************************
466 
467 !@Pawfgrtab_type
468 
469 !Retrieve sizes
470  siz_in=size(pawfgrtab_in);siz_out=size(pawfgrtab_cp)
471 
472 !Set up parallelism over atoms
473  paral_atom=(present(comm_atom));if (paral_atom) paral_atom=(xmpi_comm_size(comm_atom)>1)
474  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
475  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
476  my_atmtab_allocated=.false.
477 
478 !Determine in which case we are (parallelism, ...)
479 !No parallelism: a single copy operation
480  paral_case=0;siz_max=siz_in
481  pawfgrtab_out => pawfgrtab_cp
482  if (paral_atom) then
483    if (siz_out<siz_in) then ! Parallelism: the copy operation is a scatter
484      call get_my_natom(my_comm_atom,my_natom,siz_in)
485      if (my_natom==siz_out) then
486        call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,siz_in)
487        paral_case=1;siz_max=siz_out
488        pawfgrtab_out => pawfgrtab_cp
489      else
490        msg=' siz_out should be equal to my_natom !'
491        MSG_BUG(msg)
492      end if
493    else                            ! Parallelism: the copy operation is a gather
494      call get_my_natom(my_comm_atom,my_natom,siz_out)
495      if (my_natom==siz_in) then
496        paral_case=2;siz_max=siz_in
497      else
498        msg=' siz_in should be equal to my_natom !'
499        MSG_BUG(msg)
500      end if
501    end if
502  end if
503 
504 !First case: a simple copy or a scatter
505  if (siz_max > 0 .and. ((paral_case== 0).or.(paral_case==1))) then
506    call pawfgrtab_nullify(pawfgrtab_out)
507 
508 !  Loop on pawfgrtab components
509    do ij1=1,siz_max
510      ij=ij1; if (paral_case==1)ij=my_atmtab(ij1)
511 
512      pawfgrtab_out(ij1)%cplex=pawfgrtab_in(ij)%cplex
513      pawfgrtab_out(ij1)%expiqr_allocated=pawfgrtab_in(ij)%expiqr_allocated
514      pawfgrtab_out(ij1)%itypat=pawfgrtab_in(ij)%itypat
515      l_size=pawfgrtab_in(ij)%l_size
516      l_size2=l_size*l_size
517      pawfgrtab_out(ij1)%l_size=l_size
518      pawfgrtab_out(ij1)%gylm_allocated=pawfgrtab_in(ij)%gylm_allocated
519      pawfgrtab_out(ij1)%gylmgr_allocated=pawfgrtab_in(ij)%gylmgr_allocated
520      pawfgrtab_out(ij1)%gylmgr2_allocated=pawfgrtab_in(ij)%gylmgr2_allocated
521      pawfgrtab_out(ij1)%nfgd=pawfgrtab_in(ij)%nfgd
522      pawfgrtab_out(ij1)%nhatfr_allocated=pawfgrtab_in(ij)%nhatfr_allocated
523      pawfgrtab_out(ij1)%nhatfrgr_allocated=pawfgrtab_in(ij)%nhatfrgr_allocated
524      nspden=pawfgrtab_in(ij)%nspden
525      pawfgrtab_out(ij1)%nspden=nspden
526      pawfgrtab_out(ij1)%rfgd_allocated=pawfgrtab_in(ij)%rfgd_allocated
527      nfgd=pawfgrtab_in(ij)%nfgd
528      LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%ifftsph,(nfgd))
529      pawfgrtab_out(ij1)%ifftsph(:)=pawfgrtab_in(ij)%ifftsph(:)
530      if (pawfgrtab_out(ij1)%expiqr_allocated>0) then
531        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%expiqr,(2,nfgd))
532        pawfgrtab_out(ij1)%expiqr(:,:)=pawfgrtab_in(ij)%expiqr(:,:)
533      end if
534      if (pawfgrtab_out(ij1)%gylm_allocated>0) then
535        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylm,(nfgd,l_size2))
536        pawfgrtab_out(ij1)%gylm(:,:)=pawfgrtab_in(ij)%gylm(:,:)
537      end if
538      if (pawfgrtab_out(ij1)%gylmgr_allocated>0) then
539        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylmgr,(3,nfgd,l_size2))
540        pawfgrtab_out(ij1)%gylmgr(:,:,:)=pawfgrtab_in(ij)%gylmgr(:,:,:)
541      end if
542      if (pawfgrtab_out(ij1)%gylmgr2_allocated>0) then
543        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylmgr2,(6,nfgd,l_size2))
544        pawfgrtab_out(ij1)%gylmgr2(:,:,:)=pawfgrtab_in(ij)%gylmgr2(:,:,:)
545      end if
546      if (pawfgrtab_out(ij1)%nhatfrgr_allocated>0) then
547        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%nhatfrgr,(3,nfgd,nspden))
548        pawfgrtab_out(ij1)%nhatfrgr(:,:,:)=pawfgrtab_in(ij)%nhatfrgr(:,:,:)
549       end if
550       if (pawfgrtab_out(ij1)%rfgd_allocated>0) then
551         LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%rfgd,(3,nfgd))
552         pawfgrtab_out(ij1)%rfgd(:,:)=pawfgrtab_in(ij)%rfgd(:,:)
553       end if
554 
555     end do
556  end if
557 
558 !Second case: a gather
559  if (paral_case==2) then
560    call pawfgrtab_gather(pawfgrtab_in,pawfgrtab_cp,my_comm_atom,istat,my_atmtab)
561  end if
562 
563 !Destroy atom table used for parallelism
564  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
565 
566 end subroutine pawfgrtab_copy

m_pawfgrtab/pawfgrtab_free [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_free

FUNCTION

  Free all dynamic memory stored in a pawfgrtab datastructure

SIDE EFFECTS

  Pawfgrtab(natom) <type(pawfgrtab_type)>=atomic data given on fine rectangular grid

PARENTS

      bethe_salpeter,classify_bands,d2frnl,denfgr,exc_plot,m_fock
      m_paral_pert,m_pawfgrtab,m_wfd,pawgrnl,pawmkaewf,respfn,scfcv,screening
      sigma,wfk_analyze

CHILDREN

SOURCE

284 subroutine pawfgrtab_free(Pawfgrtab)
285 
286 
287 !This section has been created automatically by the script Abilint (TD).
288 !Do not modify the following lines by hand.
289 #undef ABI_FUNC
290 #define ABI_FUNC 'pawfgrtab_free'
291 !End of the abilint section
292 
293  implicit none
294 
295 !Arguments ------------------------------------
296 !arrays
297  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
298 
299 !Local variables-------------------------------
300 !scalars
301  integer :: iat,natom
302 
303 ! *************************************************************************
304 
305 !@Pawfgrtab_type
306 
307  natom=SIZE(Pawfgrtab);if (natom==0) return
308  do iat=1,natom
309   if (allocated(Pawfgrtab(iat)%ifftsph))  then
310     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%ifftsph)
311   end if
312   if (allocated(Pawfgrtab(iat)%gylm   ))  then
313     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylm)
314   end if
315   if (allocated(Pawfgrtab(iat)%gylmgr ))  then
316     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylmgr)
317   end if
318   if (allocated(Pawfgrtab(iat)%gylmgr2))  then
319     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylmgr2)
320   end if
321   if (allocated(Pawfgrtab(iat)%nhatfr ))  then
322     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%nhatfr)
323   end if
324   if (allocated(Pawfgrtab(iat)%nhatfrgr))  then
325     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%nhatfrgr)
326   end if
327   if (allocated(Pawfgrtab(iat)%rfgd   ))  then
328     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%rfgd)
329   end if
330   if (allocated(Pawfgrtab(iat)%expiqr))   then
331     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%expiqr)
332   end if
333   Pawfgrtab(iat)%gylm_allocated=0
334   Pawfgrtab(iat)%gylmgr_allocated=0
335   Pawfgrtab(iat)%gylmgr2_allocated=0
336   Pawfgrtab(iat)%nhatfr_allocated=0
337   Pawfgrtab(iat)%nhatfrgr_allocated=0
338   Pawfgrtab(iat)%rfgd_allocated=0
339   Pawfgrtab(iat)%expiqr_allocated=0
340  end do
341 
342 end subroutine pawfgrtab_free

m_pawfgrtab/pawfgrtab_gather [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_gather

FUNCTION

 gather a pawfgrtab
 istat =1 if pawfgrtab is yet gathered

INPUTS

  mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc
  mpicomm=MPI communicator
  pawfgrtab(:) =  pawfgrtab to be gathered

OUTPUT

  istat : 1 if yet gathered in that case, pawfgrtab_gathered contains nothing
  pawfgrtab_gathered : pawfgrtab gathered between comm_atom

PARENTS

      m_pawfgrtab,pawgrnl

CHILDREN

SOURCE

 714 subroutine pawfgrtab_gather(pawfgrtab,pawfgrtab_gathered,comm_atom,istat, &
 715 &                           mpi_atmtab) ! optional argument
 716 
 717 
 718 !This section has been created automatically by the script Abilint (TD).
 719 !Do not modify the following lines by hand.
 720 #undef ABI_FUNC
 721 #define ABI_FUNC 'pawfgrtab_gather'
 722 !End of the abilint section
 723 
 724  implicit none
 725 
 726 !Arguments ------------------------------------
 727 !scalars
 728  integer,intent(in) :: comm_atom
 729  integer, intent(out) :: istat
 730 !arrays
 731  integer,optional,target,intent(in) :: mpi_atmtab(:)
 732  type(pawfgrtab_type),intent(in) :: pawfgrtab(:)
 733  type(pawfgrtab_type),intent(inout) :: pawfgrtab_gathered(:)
 734 
 735 !Local variables-------------------------------
 736 !scalars
 737  integer :: buf_int_size,buf_int_size_all,buf_dp_size,buf_dp_size_all,i1,i2,iat,iatot
 738  integer :: ierr,ij,indx_int,indx_dp,l_size,l_size2,my_natom,natom,nfgd,nproc_atom,nspden
 739  logical :: my_atmtab_allocated,paral_atom
 740  character(len=500) :: msg
 741 !arrays
 742  integer :: bufsz(2)
 743  integer,allocatable :: buf_int(:),buf_int_all(:)
 744  integer,allocatable :: count_dp(:),count_int(:),count_tot(:),displ_dp(:),displ_int(:)
 745  integer,pointer :: my_atmtab(:)
 746  real(dp),allocatable :: buf_dp(:),buf_dp_all(:)
 747 
 748 ! *************************************************************************
 749 
 750 !@Pawfgrtab_type
 751 
 752  istat=0
 753  my_natom=size(pawfgrtab);natom=size(pawfgrtab_gathered)
 754 
 755 !Set up parallelism over atoms
 756  paral_atom=(my_natom/=natom)
 757  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
 758  my_atmtab_allocated = .False.
 759  if (paral_atom) then
 760    call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
 761  end if
 762  nproc_atom=xmpi_comm_size(comm_atom)
 763 
 764 !Without parallelism, just copy input to output
 765  if (.not.paral_atom) then
 766    istat=1;return
 767  end if
 768 
 769 !Compute size of buffers
 770  buf_int_size=0
 771  do iat=1,my_natom
 772    buf_int_size=buf_int_size+pawfgrtab(iat)%nfgd+13
 773  end do
 774  buf_dp_size=0
 775  do iat=1,my_natom
 776    if (pawfgrtab(iat)%expiqr_allocated>0) then
 777      buf_dp_size=buf_dp_size + size(pawfgrtab(iat)%expiqr(:,:))
 778    end if
 779    if (pawfgrtab(iat)%gylm_allocated>0) then
 780      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylm(:,:))
 781    end if
 782    if (pawfgrtab(iat)%gylmgr_allocated>0) then
 783      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylmgr(:,:,:))
 784    end if
 785    if (pawfgrtab(iat)%gylmgr2_allocated>0) then
 786      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylmgr2(:,:,:))
 787    end if
 788    if (pawfgrtab(iat)%nhatfr_allocated>0) then
 789      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%nhatfr(:,:))
 790    end if
 791    if (pawfgrtab(iat)%nhatfrgr_allocated>0) then
 792      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%nhatfrgr(:,:,:))
 793    end if
 794    if (pawfgrtab(iat)%rfgd_allocated>0) then
 795      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%rfgd(:,:))
 796    end if
 797  end do
 798 
 799 !Fill in input buffers
 800  LIBPAW_ALLOCATE(buf_int,(buf_int_size))
 801  LIBPAW_ALLOCATE(buf_dp,(buf_dp_size))
 802  indx_int=1;indx_dp=1
 803  do iat=1,my_natom
 804    l_size=pawfgrtab(iat)%l_size
 805    nfgd=pawfgrtab(iat)%nfgd
 806    nspden=pawfgrtab(iat)%nspden
 807    buf_int(indx_int)=my_atmtab(iat); indx_int=indx_int+1;
 808    buf_int(indx_int)=pawfgrtab(iat)%cplex; indx_int=indx_int+1;
 809    buf_int(indx_int)=pawfgrtab(iat)%itypat; indx_int=indx_int+1;
 810    buf_int(indx_int)=pawfgrtab(iat)%expiqr_allocated; indx_int=indx_int+1;
 811    buf_int(indx_int)=pawfgrtab(iat)%l_size; indx_int=indx_int+1;
 812    buf_int(indx_int)=pawfgrtab(iat)%gylm_allocated; indx_int=indx_int+1;
 813    buf_int(indx_int)=pawfgrtab(iat)%gylmgr_allocated; indx_int=indx_int+1;
 814    buf_int(indx_int)=pawfgrtab(iat)%gylmgr2_allocated; indx_int=indx_int+1;
 815    buf_int(indx_int)=pawfgrtab(iat)%nfgd; indx_int=indx_int+1;
 816    buf_int(indx_int)=pawfgrtab(iat)%nhatfr_allocated; indx_int=indx_int+1;
 817    buf_int(indx_int)=pawfgrtab(iat)%nhatfrgr_allocated; indx_int=indx_int+1;
 818    buf_int(indx_int)=pawfgrtab(iat)%nspden; indx_int=indx_int+1;
 819    buf_int(indx_int)=pawfgrtab(iat)%rfgd_allocated; indx_int=indx_int+1;
 820    if (nfgd>0) then
 821      buf_int(indx_int:indx_int+nfgd-1)=pawfgrtab(iat)%ifftsph(1:nfgd)
 822      indx_int=indx_int+nfgd;
 823    end if
 824    if (pawfgrtab(iat)%expiqr_allocated>0) then
 825      do i1=1,nfgd
 826        buf_dp(indx_dp:indx_dp+1)=pawfgrtab(iat)%expiqr(1:2,i1)
 827        indx_dp=indx_dp+2
 828      end do
 829    end if
 830    l_size2=l_size*l_size
 831    if (pawfgrtab(iat)%gylm_allocated>0) then
 832      do i1=1,l_size2
 833        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab(iat)%gylm(1:nfgd,i1)
 834       indx_dp=indx_dp+nfgd
 835      end do
 836    end if
 837    if (pawfgrtab(iat)%gylmgr_allocated>0) then
 838      do i1=1,l_size2
 839        do i2=1,nfgd
 840          buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%gylmgr(1:3,i2,i1)
 841          indx_dp=indx_dp+3
 842        end do
 843      end do
 844    end if
 845    if (pawfgrtab(iat)%gylmgr2_allocated>0) then
 846      do i1=1,l_size2
 847        do i2=1,nfgd
 848          buf_dp(indx_dp:indx_dp+5)=pawfgrtab(iat)%gylmgr2(1:6,i2,i1)
 849          indx_dp=indx_dp+6
 850        end do
 851      end do
 852    end if
 853    if (pawfgrtab(iat)%nhatfr_allocated>0) then
 854      do i1=1,nspden
 855        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab(iat)%nhatfr(1:nfgd,i1)
 856        indx_dp=indx_dp+nfgd
 857      end do
 858    end if
 859    if (pawfgrtab(iat)%nhatfrgr_allocated>0) then
 860      do i1=1,nspden
 861        do i2=1,nfgd
 862          buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%nhatfrgr(1:3,i2,i1)
 863          indx_dp=indx_dp+3
 864        end do
 865      end do
 866    end if
 867    if (pawfgrtab(iat)%rfgd_allocated>0) then
 868      do i1=1,nfgd
 869        buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%rfgd(1:3,i1)
 870        indx_dp=indx_dp+3
 871      end do
 872    end if
 873  end do
 874  if (indx_int/=1+buf_int_size) then
 875    msg='Error (1) in pawfgrtab_gather: wrong buffer sizes !'
 876    MSG_BUG(msg)
 877  end if
 878  if (indx_dp/=1+buf_dp_size) then
 879    msg='Error (2) in pawfgrtab_gather: wrong buffer sizes !'
 880    MSG_BUG(msg)
 881  end if
 882 
 883 !Communicate (1 gather for integers, 1 gather for reals)
 884  LIBPAW_ALLOCATE(count_int,(nproc_atom))
 885  LIBPAW_ALLOCATE(displ_int,(nproc_atom))
 886  LIBPAW_ALLOCATE(count_dp ,(nproc_atom))
 887  LIBPAW_ALLOCATE(displ_dp ,(nproc_atom))
 888  LIBPAW_ALLOCATE(count_tot,(2*nproc_atom))
 889  bufsz(1)=buf_int_size; bufsz(2)=buf_dp_size
 890  call xmpi_allgather(bufsz,2,count_tot,comm_atom,ierr)
 891  do ij=1,nproc_atom
 892    count_int(ij)=count_tot(2*ij-1)
 893    count_dp (ij)=count_tot(2*ij)
 894  end do
 895  displ_int(1)=0;displ_dp(1)=0
 896  do ij=2,nproc_atom
 897    displ_int(ij)=displ_int(ij-1)+count_int(ij-1)
 898    displ_dp (ij)=displ_dp (ij-1)+count_dp (ij-1)
 899  end do
 900  buf_int_size_all=sum(count_int)
 901  buf_dp_size_all =sum(count_dp)
 902  LIBPAW_DEALLOCATE(count_tot)
 903  LIBPAW_ALLOCATE(buf_int_all,(buf_int_size_all))
 904  LIBPAW_ALLOCATE(buf_dp_all ,(buf_dp_size_all))
 905  call xmpi_allgatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,comm_atom,ierr)
 906  call xmpi_allgatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,comm_atom,ierr)
 907  LIBPAW_DEALLOCATE(count_int)
 908  LIBPAW_DEALLOCATE(displ_int)
 909  LIBPAW_DEALLOCATE(count_dp)
 910  LIBPAW_DEALLOCATE(displ_dp)
 911 
 912 !Retrieve output datastructures
 913  indx_int=1; indx_dp=1
 914  call pawfgrtab_free(pawfgrtab_gathered)
 915  call pawfgrtab_nullify(pawfgrtab_gathered)
 916  do iat=1,natom
 917    iatot=buf_int_all(indx_int); indx_int=indx_int+1
 918    pawfgrtab_gathered(iatot)%cplex=buf_int_all(indx_int); indx_int=indx_int+1;
 919    pawfgrtab_gathered(iatot)%itypat=buf_int_all(indx_int); indx_int=indx_int+1;
 920    pawfgrtab_gathered(iatot)%expiqr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 921    pawfgrtab_gathered(iatot)%l_size=buf_int_all(indx_int); indx_int=indx_int+1;
 922    pawfgrtab_gathered(iatot)%gylm_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 923    pawfgrtab_gathered(iatot)%gylmgr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 924    pawfgrtab_gathered(iatot)%gylmgr2_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 925    pawfgrtab_gathered(iatot)%nfgd=buf_int_all(indx_int); indx_int=indx_int+1;
 926    pawfgrtab_gathered(iatot)%nhatfr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 927    pawfgrtab_gathered(iatot)%nhatfrgr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 928    pawfgrtab_gathered(iatot)%nspden=buf_int_all(indx_int); indx_int=indx_int+1;
 929    pawfgrtab_gathered(iatot)%rfgd_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
 930    nfgd=pawfgrtab_gathered(iatot)%nfgd
 931    l_size=pawfgrtab_gathered(iatot)%l_size
 932    nspden= pawfgrtab_gathered(iatot)%nspden
 933    l_size2=l_size*l_size
 934    if (nfgd>0) then
 935      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%ifftsph,(nfgd))
 936      pawfgrtab_gathered(iatot)%ifftsph(1:nfgd)=buf_int_all(indx_int:indx_int+nfgd-1)
 937      indx_int=indx_int+nfgd
 938    end if
 939    if (pawfgrtab_gathered(iatot)%expiqr_allocated>0) then
 940      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%expiqr,(2,nfgd))
 941      do i1=1,nfgd
 942        pawfgrtab_gathered(iatot)%expiqr(1:2,i1)= buf_dp_all(indx_dp:indx_dp+1)
 943        indx_dp=indx_dp+2
 944      end do
 945    end if
 946    if (pawfgrtab_gathered(iatot)%gylm_allocated>0) then
 947      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylm,(nfgd,l_size2))
 948      do i1=1,l_size2
 949        pawfgrtab_gathered(iatot)%gylm(1:nfgd,i1)=buf_dp_all(indx_dp:indx_dp+nfgd-1)
 950        indx_dp=indx_dp+nfgd
 951      end do
 952    end if
 953    if (pawfgrtab_gathered(iatot)%gylmgr_allocated>0) then
 954      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylmgr,(3,nfgd,l_size2))
 955      do i1=1,l_size2
 956        do i2=1,nfgd
 957          pawfgrtab_gathered(iatot)%gylmgr(1:3,i2,i1)=buf_dp_all(indx_dp:indx_dp+2)
 958          indx_dp=indx_dp+3
 959        end do
 960      end do
 961    end if
 962    if (pawfgrtab_gathered(iatot)%gylmgr2_allocated>0) then
 963      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylmgr2,(6,nfgd,l_size2))
 964      do i1=1,l_size2
 965        do i2=1,nfgd
 966          pawfgrtab_gathered(iatot)%gylmgr2(1:6,i2,i1)=buf_dp_all(indx_dp:indx_dp+5)
 967          indx_dp=indx_dp+6
 968        end do
 969      end do
 970    end if
 971    if (pawfgrtab_gathered(iatot)%nhatfr_allocated>0) then
 972      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%nhatfr,(nfgd,nspden))
 973      do i1=1,nspden
 974        pawfgrtab_gathered(iatot)%nhatfr(1:nfgd,i1)=buf_dp_all(indx_dp:indx_dp+nfgd-1)
 975        indx_dp=indx_dp+nfgd
 976      end do
 977    end if
 978    if (pawfgrtab_gathered(iatot)%nhatfrgr_allocated>0) then
 979      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%nhatfrgr,(3,nfgd,nspden))
 980      do i1=1,nspden
 981        do i2=1,nfgd
 982          pawfgrtab_gathered(iatot)%nhatfrgr(1:3,i2,i1)=buf_dp_all(indx_dp:indx_dp+2)
 983          indx_dp=indx_dp+3
 984        end do
 985      end do
 986    end if
 987    if (pawfgrtab_gathered(iatot)%rfgd_allocated>0) then
 988      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%rfgd,(3,nfgd))
 989      do i1=1,nfgd
 990        pawfgrtab_gathered(iatot)%rfgd(1:3,i1)=buf_dp_all(indx_dp:indx_dp+2)
 991        indx_dp=indx_dp+3
 992      end do
 993    end if
 994  end do
 995 
 996 !Free memory
 997  LIBPAW_DEALLOCATE(buf_int)
 998  LIBPAW_DEALLOCATE(buf_int_all)
 999  LIBPAW_DEALLOCATE(buf_dp)
1000  LIBPAW_DEALLOCATE(buf_dp_all)
1001 
1002 !Destroy atom table
1003  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
1004 
1005 end subroutine pawfgrtab_gather

m_pawfgrtab/pawfgrtab_init [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_init

FUNCTION

  Initialize a pawfgrtab datastructure

OUTPUT

  Pawfgrtab(natom) <type(pawfgrtab_type)>=atomic data given on fine rectangular grid

PARENTS

      bethe_salpeter,classify_bands,d2frnl,denfgr,exc_plot,m_fock,m_wfd
      pawmkaewf,respfn,scfcv,screening,sigma,wfk_analyze

CHILDREN

SOURCE

184 subroutine pawfgrtab_init(Pawfgrtab,cplex,l_size_atm,nspden,typat,&
185 &                         mpi_atmtab,comm_atom) ! optional arguments (parallelism)
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 'pawfgrtab_init'
192 !End of the abilint section
193 
194  implicit none
195 
196 !Arguments ------------------------------------
197 !scalars
198  integer,intent(in) :: cplex,nspden
199 !arrays
200  integer,intent(in) :: l_size_atm(:),typat(:)
201  integer,optional,intent(in) :: comm_atom
202  integer,optional,target,intent(in) :: mpi_atmtab(:)
203  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
204 
205 !Local variables-------------------------------
206 !scalars
207  integer :: iat,iatom_tot,my_comm_atom,my_natom,natom
208  logical :: my_atmtab_allocated,paral_atom
209  character(len=500) :: msg
210 !arrays
211  integer,pointer :: my_atmtab(:)
212 
213 ! *************************************************************************
214 
215 !@Pawfgrtab_type
216 
217  my_natom=SIZE(Pawfgrtab);if (my_natom==0) return
218  natom=SIZE(typat)
219  if (my_natom/=SIZE(l_size_atm)) then
220   msg='Sizes of assumed shape arrays do not match'
221   MSG_BUG(msg)
222  end if
223 
224 !Set up parallelism over atoms
225  paral_atom=(present(comm_atom).and.(my_natom/=natom))
226  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
227  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
228  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
229 
230  call pawfgrtab_nullify(Pawfgrtab)
231 
232  do iat=1,my_natom
233   iatom_tot=iat; if (paral_atom) iatom_tot=my_atmtab(iat)
234 
235   Pawfgrtab(iat)%cplex             = cplex
236   Pawfgrtab(iat)%itypat            = typat(iatom_tot)
237   Pawfgrtab(iat)%nspden            = nspden
238   Pawfgrtab(iat)%l_size            = l_size_atm(iat)
239   Pawfgrtab(iat)%nfgd              = 0
240   LIBPAW_ALLOCATE(Pawfgrtab(iat)%ifftsph,(0))
241   Pawfgrtab(iat)%gylm_allocated    = 0
242   LIBPAW_ALLOCATE(Pawfgrtab(iat)%gylm,(0,0))
243   Pawfgrtab(iat)%gylmgr_allocated  = 0
244   LIBPAW_ALLOCATE(Pawfgrtab(iat)%gylmgr,(0,0,0))
245   Pawfgrtab(iat)%gylmgr2_allocated = 0
246   LIBPAW_ALLOCATE(Pawfgrtab(iat)%gylmgr2,(0,0,0))
247   Pawfgrtab(iat)%nhatfr_allocated  = 0
248   LIBPAW_ALLOCATE(Pawfgrtab(iat)%nhatfr,(0,0))
249   Pawfgrtab(iat)%nhatfrgr_allocated  = 0
250   LIBPAW_ALLOCATE(Pawfgrtab(iat)%nhatfrgr,(0,0,0))
251   Pawfgrtab(iat)%rfgd_allocated    = 0
252   LIBPAW_ALLOCATE(Pawfgrtab(iat)%rfgd,(0,0))
253   Pawfgrtab(iat)%expiqr_allocated  = 0
254   LIBPAW_ALLOCATE(Pawfgrtab(iat)%expiqr,(0,0))
255  end do
256 
257 !Destroy atom table used for parallelism
258  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
259 
260 end subroutine pawfgrtab_init

m_pawfgrtab/pawfgrtab_isendreceive_fillbuffer [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

  pawfgrtab_isendreceive_fillbuffer

FUNCTION

  Extract from pawfgrtab and from the global index of atoms
  the buffers to send in a sending operation
  This function has to be coupled with a call to pawfgrtab_isendreceive_getbuffer

INPUTS

  atm_indx_send(1:total number of atoms)= array for send operation,
                 Given an index of atom in global numbering, give its index
                 in the table of atoms treated by current processor
                 or -1 if the atoms is not treated by current processor
  npawfgrtab_send= number of sent atoms
  pawfgrtab= data structure from which are extract buffer int and buffer dp

OUTPUT

  buf_int= buffer of integers to be sent
  buf_int_size= size of buffer of integers
  buf_dp= buffer of double precision numbers to be sent
  buf_dp_size= size of buffer of double precision numbers

PARENTS

      m_pawfgrtab

CHILDREN

SOURCE

1552 subroutine pawfgrtab_isendreceive_fillbuffer(pawfgrtab, atmtab_send,atm_indx_send,npawfgrtab_send,&
1553 &                                            buf_int,buf_int_size,buf_dp,buf_dp_size)
1554 
1555 
1556 !This section has been created automatically by the script Abilint (TD).
1557 !Do not modify the following lines by hand.
1558 #undef ABI_FUNC
1559 #define ABI_FUNC 'pawfgrtab_isendreceive_fillbuffer'
1560 !End of the abilint section
1561 
1562 implicit none
1563 
1564 !Arguments ------------------------------------
1565 !scalars
1566  integer,intent(out) :: buf_int_size,buf_dp_size
1567  integer,intent(in) :: npawfgrtab_send
1568 !arrays
1569  integer,intent(in) :: atmtab_send(:),atm_indx_send(:)
1570  integer,intent(out),allocatable :: buf_int(:)
1571  real(dp),intent(out),allocatable :: buf_dp(:)
1572  type(pawfgrtab_type),target,intent(inout) :: pawfgrtab(:)
1573 
1574 !Local variables-------------------------------
1575 !scalars
1576  integer :: i1,i2,iat,iat1,iatom_tot,indx_int,indx_dp,l_size,l_size2,nfgd,nspden
1577  character(len=500) :: msg
1578  type(pawfgrtab_type),pointer :: pawfgrtab1
1579 
1580 ! *********************************************************************
1581 
1582 !Compute size of buffers
1583  buf_int_size=0;buf_dp_size=0
1584  do iat1=1,npawfgrtab_send
1585    iatom_tot=atmtab_send(iat1)
1586    iat=atm_indx_send(iatom_tot)
1587    pawfgrtab1=>pawfgrtab(iat)
1588    buf_int_size=buf_int_size+pawfgrtab1%nfgd+13
1589    if (pawfgrtab1%expiqr_allocated>0) then
1590      buf_dp_size=buf_dp_size + size(pawfgrtab1%expiqr(:,:))
1591    end if
1592    if (pawfgrtab1%gylm_allocated>0) then
1593      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylm(:,:))
1594    end if
1595    if (pawfgrtab1%gylmgr_allocated>0) then
1596      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylmgr(:,:,:))
1597    end if
1598    if (pawfgrtab1%gylmgr2_allocated>0 ) then
1599      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylmgr2(:,:,:))
1600    end if
1601    if (pawfgrtab1%nhatfr_allocated>0) then
1602      buf_dp_size=buf_dp_size+size(pawfgrtab1%nhatfr(:,:))
1603    end if
1604    if (pawfgrtab1%nhatfrgr_allocated>0) then
1605      buf_dp_size=buf_dp_size+size(pawfgrtab1%nhatfrgr(:,:,:))
1606    end if
1607    if (pawfgrtab1%rfgd_allocated>0) then
1608      buf_dp_size=buf_dp_size+size(pawfgrtab1%rfgd(:,:))
1609    end if
1610  end do
1611 
1612 !Fill in input buffers
1613  LIBPAW_ALLOCATE(buf_int,(buf_int_size))
1614  LIBPAW_ALLOCATE(buf_dp,(buf_dp_size))
1615  indx_int=1;indx_dp=1
1616  do iat1=1,npawfgrtab_send
1617     iatom_tot=atmtab_send(iat1)
1618    iat=atm_indx_send(iatom_tot)
1619    pawfgrtab1=>pawfgrtab(iat)
1620    l_size=pawfgrtab1%l_size
1621    nfgd=pawfgrtab1%nfgd
1622    nspden=pawfgrtab1%nspden
1623    buf_int(indx_int)=iatom_tot; indx_int=indx_int+1;
1624    buf_int(indx_int)=pawfgrtab1%cplex; indx_int=indx_int+1;
1625    buf_int(indx_int)=pawfgrtab1%itypat; indx_int=indx_int+1;
1626    buf_int(indx_int)=pawfgrtab1%expiqr_allocated; indx_int=indx_int+1;
1627    buf_int(indx_int)=pawfgrtab1%l_size; indx_int=indx_int+1;
1628    buf_int(indx_int)=pawfgrtab1%gylm_allocated; indx_int=indx_int+1;
1629    buf_int(indx_int)=pawfgrtab1%gylmgr_allocated; indx_int=indx_int+1;
1630    buf_int(indx_int)=pawfgrtab1%gylmgr2_allocated; indx_int=indx_int+1;
1631    buf_int(indx_int)=pawfgrtab1%nfgd; indx_int=indx_int+1;
1632    buf_int(indx_int)=pawfgrtab1%nhatfr_allocated; indx_int=indx_int+1;
1633    buf_int(indx_int)=pawfgrtab1%nhatfrgr_allocated; indx_int=indx_int+1;
1634    buf_int(indx_int)=pawfgrtab1%nspden; indx_int=indx_int+1;
1635    buf_int(indx_int)=pawfgrtab1%rfgd_allocated; indx_int=indx_int+1;
1636    if (nfgd>0) then
1637      buf_int(indx_int:indx_int+nfgd-1)=pawfgrtab1%ifftsph(1:nfgd)
1638      indx_int=indx_int+nfgd;
1639    end if
1640    if (pawfgrtab1%expiqr_allocated>0) then
1641      do i1=1,nfgd
1642        buf_dp(indx_dp:indx_dp+1)=pawfgrtab1%expiqr(1:2,i1)
1643        indx_dp=indx_dp+2
1644      end do
1645    end if
1646    l_size2=l_size*l_size
1647    if (pawfgrtab1%gylm_allocated>0) then
1648      do i1=1,l_size2
1649        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab1%gylm(1:nfgd,i1)
1650       indx_dp=indx_dp+nfgd
1651      end do
1652    end if
1653    if (pawfgrtab1%gylmgr_allocated>0) then
1654      do i1=1,l_size2
1655        do i2=1,nfgd
1656          buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%gylmgr(1:3,i2,i1)
1657          indx_dp=indx_dp+3
1658        end do
1659      end do
1660    end if
1661    if (pawfgrtab1%gylmgr2_allocated>0) then
1662      do i1=1,l_size2
1663        do i2=1,nfgd
1664          buf_dp(indx_dp:indx_dp+5)=pawfgrtab1%gylmgr2(1:6,i2,i1)
1665          indx_dp=indx_dp+6
1666        end do
1667      end do
1668    end if
1669    if (pawfgrtab1%nhatfr_allocated>0) then
1670      do i1=1,nspden
1671        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab1%nhatfr(1:nfgd,i1)
1672        indx_dp=indx_dp+nfgd
1673      end do
1674    end if
1675    if (pawfgrtab1%nhatfrgr_allocated>0) then
1676      do i1=1,nspden
1677        do i2=1,nfgd
1678          buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%nhatfrgr(1:3,i2,i1)
1679          indx_dp=indx_dp+3
1680        end do
1681      end do
1682    end if
1683    if (pawfgrtab1%rfgd_allocated>0) then
1684      do i1=1,nfgd
1685        buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%rfgd(1:3,i1)
1686        indx_dp=indx_dp+3
1687      end do
1688    end if
1689  end do
1690  if ((indx_int-1/=buf_int_size).or.(indx_dp-1/=buf_dp_size)) then
1691    write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size
1692    MSG_BUG(msg)
1693  end if
1694 
1695 end subroutine pawfgrtab_isendreceive_fillbuffer

m_pawfgrtab/pawfgrtab_isendreceive_getbuffer [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

  pawfgrtab_isendreceive_getbuffer

FUNCTION

  Fill a pawfgrtab structure with the buffers received in a receive operation
  This buffer should have been first extracted by a call to pawfgrtab_isendreceive_fillbuffer

INPUTS

  atm_indx_recv(1:total number of atoms)= array for receive operation
                 Given an index of atom in global numbering, give its index
                 in the table of atoms treated by current processor
                 or -1 if the atoms is not treated by current processor
  buf_int= buffer of receive integers
  buf_dp= buffer of receive double precision numbers
  npawfgrtab_send= number of sent atoms

OUTPUT

  pawfgrtab= output datastructure filled with buffers receive in a receive operation

PARENTS

      m_pawfgrtab

CHILDREN

SOURCE

1400 subroutine pawfgrtab_isendreceive_getbuffer(pawfgrtab,npawfgrtab_send,atm_indx_recv,buf_int,buf_dp)
1401 
1402 
1403 !This section has been created automatically by the script Abilint (TD).
1404 !Do not modify the following lines by hand.
1405 #undef ABI_FUNC
1406 #define ABI_FUNC 'pawfgrtab_isendreceive_getbuffer'
1407 !End of the abilint section
1408 
1409 implicit none
1410 
1411 !Arguments ------------------------------------
1412 !scalars
1413  integer,intent(in) :: npawfgrtab_send
1414 !arrays
1415  integer,intent(in) ::atm_indx_recv(:),buf_int(:)
1416  real(dp),intent(in) :: buf_dp(:)
1417  type(pawfgrtab_type),target,intent(inout) :: pawfgrtab(:)
1418 
1419 !Local variables-------------------------------
1420 !scalars
1421  integer :: buf_dp_size,buf_int_size,i1,i2,iat,iatom_tot,ij,indx_int,indx_dp
1422  integer :: l_size,l_size2,nfgd,nspden
1423  character(len=500) :: msg
1424  type(pawfgrtab_type),pointer :: pawfgrtab1
1425 
1426 ! *********************************************************************
1427 
1428  buf_int_size=size(buf_int)
1429  buf_dp_size=size(buf_dp)
1430  indx_int=1; indx_dp=1
1431 
1432  do ij=1,npawfgrtab_send
1433    iatom_tot=buf_int(indx_int); indx_int=indx_int+1
1434    iat= atm_indx_recv(iatom_tot)
1435    pawfgrtab1=>pawfgrtab(iat)
1436    pawfgrtab1%cplex=buf_int(indx_int); indx_int=indx_int+1;
1437    pawfgrtab1%itypat=buf_int(indx_int); indx_int=indx_int+1;
1438    pawfgrtab1%expiqr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1439    pawfgrtab1%l_size=buf_int(indx_int); indx_int=indx_int+1;
1440    pawfgrtab1%gylm_allocated=buf_int(indx_int); indx_int=indx_int+1;
1441    pawfgrtab1%gylmgr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1442    pawfgrtab1%gylmgr2_allocated=buf_int(indx_int); indx_int=indx_int+1;
1443    pawfgrtab1%nfgd=buf_int(indx_int); indx_int=indx_int+1;
1444    pawfgrtab1%nhatfr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1445    pawfgrtab1%nhatfrgr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1446    pawfgrtab1%nspden=buf_int(indx_int); indx_int=indx_int+1;
1447    pawfgrtab1%rfgd_allocated=buf_int(indx_int); indx_int=indx_int+1;
1448    nfgd=pawfgrtab1%nfgd
1449    l_size=pawfgrtab1%l_size
1450    nspden= pawfgrtab1%nspden
1451    l_size2=l_size*l_size
1452    if (nfgd>0) then
1453      LIBPAW_ALLOCATE(pawfgrtab1%ifftsph,(nfgd))
1454      pawfgrtab1%ifftsph(1:nfgd)=buf_int(indx_int:indx_int+nfgd-1)
1455      indx_int=indx_int+nfgd
1456    end if
1457    if (pawfgrtab1%expiqr_allocated>0) then
1458      LIBPAW_ALLOCATE(pawfgrtab1%expiqr,(2,nfgd))
1459      do i1=1,nfgd
1460        pawfgrtab1%expiqr(1:2,i1)= buf_dp(indx_dp:indx_dp+1)
1461        indx_dp=indx_dp+2
1462      end do
1463    end if
1464    if (pawfgrtab1%gylm_allocated>0) then
1465      LIBPAW_ALLOCATE(pawfgrtab1%gylm,(nfgd,l_size2))
1466      do i1=1,l_size2
1467        pawfgrtab1%gylm(1:nfgd,i1)=buf_dp(indx_dp:indx_dp+nfgd-1)
1468        indx_dp=indx_dp+nfgd
1469      end do
1470    end if
1471    if (pawfgrtab1%gylmgr_allocated>0) then
1472      LIBPAW_ALLOCATE(pawfgrtab1%gylmgr,(3,nfgd,l_size2))
1473      do i1=1,l_size2
1474        do i2=1,nfgd
1475          pawfgrtab1%gylmgr(1:3,i2,i1)=buf_dp(indx_dp:indx_dp+2)
1476          indx_dp=indx_dp+3
1477        end do
1478      end do
1479    end if
1480    if (pawfgrtab1%gylmgr2_allocated>0) then
1481      LIBPAW_ALLOCATE(pawfgrtab1%gylmgr2,(6,nfgd,l_size2))
1482      do i1=1,l_size2
1483        do i2=1,nfgd
1484          pawfgrtab1%gylmgr2(1:6,i2,i1)=buf_dp(indx_dp:indx_dp+5)
1485          indx_dp=indx_dp+6
1486        end do
1487      end do
1488    end if
1489    if (pawfgrtab1%nhatfr_allocated>0) then
1490      LIBPAW_ALLOCATE(pawfgrtab1%nhatfr,(nfgd,nspden))
1491      do i1=1,nspden
1492        pawfgrtab1%nhatfr(1:nfgd,i1)=buf_dp(indx_dp:indx_dp+nfgd-1)
1493        indx_dp=indx_dp+nfgd
1494      end do
1495    end if
1496    if (pawfgrtab1%nhatfrgr_allocated>0) then
1497      LIBPAW_ALLOCATE(pawfgrtab1%nhatfrgr,(3,nfgd,nspden))
1498      do i1=1,nspden
1499        do i2=1,nfgd
1500          pawfgrtab1%nhatfrgr(1:3,i2,i1)=buf_dp(indx_dp:indx_dp+2)
1501          indx_dp=indx_dp+3
1502        end do
1503      end do
1504    end if
1505    if (pawfgrtab1%rfgd_allocated>0) then
1506      LIBPAW_ALLOCATE(pawfgrtab1%rfgd,(3,nfgd))
1507      do i1=1,nfgd
1508        pawfgrtab1%rfgd(1:3,i1)=buf_dp(indx_dp:indx_dp+2)
1509        indx_dp=indx_dp+3
1510      end do
1511    end if
1512  end do
1513  if ((indx_int/=1+buf_int_size).or.(indx_dp/=1+buf_dp_size)) then
1514    write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size
1515    MSG_BUG(msg)
1516  end if
1517 
1518 end subroutine pawfgrtab_isendreceive_getbuffer

m_pawfgrtab/pawfgrtab_nullify [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_nullify

FUNCTION

  Nullify the pointers in a pawfgrtab datastructure

SIDE EFFECTS

  Pawfgrtab(:) <type(pawfgrtab_type)>=atomic data given on fine rectangular grid

PARENTS

      m_fock,m_pawfgrtab,pawgrnl

CHILDREN

SOURCE

364 subroutine pawfgrtab_nullify(Pawfgrtab)
365 
366 
367 !This section has been created automatically by the script Abilint (TD).
368 !Do not modify the following lines by hand.
369 #undef ABI_FUNC
370 #define ABI_FUNC 'pawfgrtab_nullify'
371 !End of the abilint section
372 
373  implicit none
374 
375 !Arguments ------------------------------------
376 !arrays
377  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
378 
379 !Local variables-------------------------------
380 !scalars
381  integer :: iat,natom
382 
383 ! *************************************************************************
384 
385 !@Pawfgrtab_type
386 
387 ! MGPAW: This one could be removed/renamed,
388 ! variables can be initialized in the datatype declaration
389 ! Do we need to expose this in the public API?
390 
391  natom=SIZE(Pawfgrtab);if (natom==0) return
392  do iat=1,natom
393   Pawfgrtab(iat)%nfgd              = 0
394   Pawfgrtab(iat)%gylm_allocated    = 0
395   Pawfgrtab(iat)%gylmgr_allocated  = 0
396   Pawfgrtab(iat)%gylmgr2_allocated = 0
397   Pawfgrtab(iat)%nhatfr_allocated  = 0
398   Pawfgrtab(iat)%nhatfrgr_allocated= 0
399   Pawfgrtab(iat)%rfgd_allocated    = 0
400   Pawfgrtab(iat)%expiqr_allocated  = 0
401  end do
402 
403 end subroutine pawfgrtab_nullify

m_pawfgrtab/pawfgrtab_print [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_print

FUNCTION

  Reports basic info on the pawfgrtab datatype.

INPUTS

  Pawfgrtab<pawfgrtab_type>=The datatype to be printed
  [mode_paral]=either "COLL" or "PERS", "COLL" if None.
  [unit]=Unit number for output, std_out if None.
  [prtvol]=Verbosity level, lowest if None.
  [mpi_atmtab(:)]=indexes of the atoms treated by current proc
  [comm_atom]=MPI communicator over atoms
  [natom]=total number of atom (needed if parallelism over atoms is activated)
          if Pawfgrtab is distributed, natom is different from size(Pawfgrtab).

OUTPUT

 (only writing)

PARENTS

      exc_plot,m_wfd,pawmkaewf,sigma,wfk_analyze

CHILDREN

SOURCE

598 subroutine pawfgrtab_print(Pawfgrtab,natom,unit,prtvol,mode_paral,mpi_atmtab,comm_atom)
599 
600 
601 !This section has been created automatically by the script Abilint (TD).
602 !Do not modify the following lines by hand.
603 #undef ABI_FUNC
604 #define ABI_FUNC 'pawfgrtab_print'
605 !End of the abilint section
606 
607  implicit none
608 
609 !Arguments ------------------------------------
610 !scalars
611  integer,optional,intent(in) :: comm_atom,natom,prtvol,unit
612  character(len=4),optional,intent(in) :: mode_paral
613 !arrays
614  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
615  integer,optional,target,intent(in)  :: mpi_atmtab(:)
616 
617 !Local variables-------------------------------
618 !scalars
619  integer :: iat,iatom_tot,my_comm_atom,my_natom,my_unt,my_prtvol,size_Pawfgrtab
620  logical :: my_atmtab_allocated,paral_atom
621  character(len=4) :: my_mode
622  character(len=500) :: msg
623 !arrays
624  integer,pointer :: my_atmtab(:)
625 
626 ! *************************************************************************
627 
628 !@Pawfgrtab_type
629 
630  size_Pawfgrtab=SIZE(Pawfgrtab);if (size_Pawfgrtab==0) return
631 
632  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
633  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
634  my_mode  ='PERS' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
635  my_natom=size_Pawfgrtab; if (PRESENT(natom))  my_natom=natom
636 
637 !Set up parallelism over atoms
638  paral_atom=(present(comm_atom).and.my_natom/=size_Pawfgrtab)
639  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
640  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
641  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,&
642 &                   my_natom,my_natom_ref=size_Pawfgrtab)
643 
644  write(msg,'(3a)')ch10,' === Content of the pawfgrtab datatype === ',ch10
645  call wrtout(my_unt,msg,my_mode)
646  do iat=1,my_natom
647    iatom_tot=iat;if(paral_atom) iatom_tot=my_atmtab(iat)
648 
649    write(msg,'(6(2a,i0))')ch10,&
650 &    ' > For atom number : ',iatom_tot,ch10,&
651 &    '    cplex= ',Pawfgrtab(iat)%cplex,ch10,&
652 &    '    # of spins components for density= ',Pawfgrtab(iat)%nspden,ch10,&
653 &    '    Type of atom= ',Pawfgrtab(iat)%itypat,ch10,&
654 &    '    1+ Max l in Gaunt coefficients= ',Pawfgrtab(iat)%l_size,ch10,&
655 &    '    Number of fine FFT points in PAW sphere= ',Pawfgrtab(iat)%nfgd
656    call wrtout(my_unt,msg,my_mode)
657 
658    if (my_prtvol>=3) then
659      write(msg,'(a,7(a,i2,a))')ch10,&
660 &      '    rfgd_allocated    : ',Pawfgrtab(iat)%rfgd_allocated,ch10,&
661 &      '    gylm_allocated    : ',Pawfgrtab(iat)%gylm_allocated,ch10,&
662 &      '    gylmgr_allocated  : ',Pawfgrtab(iat)%gylmgr_allocated,ch10,&
663 &      '    gylmgr2_allocated : ',Pawfgrtab(iat)%gylmgr2_allocated,ch10,&
664 &      '    nhatgr_allocated  : ',Pawfgrtab(iat)%nhatfr_allocated,ch10,&
665 &      '    nhatfrgr_allocated: ',Pawfgrtab(iat)%nhatfrgr_allocated,ch10,&
666 &      '    expiqr_allocated  : ',Pawfgrtab(iat)%expiqr_allocated,ch10
667      call wrtout(my_unt,msg,my_mode)
668 
669    end if
670 
671 !  These huge arrays are not printed out!
672 !  Pawfgrtab(iat)%ifftsph
673 !  Pawfgrtab(iat)%rfgd
674 !  Pawfgrtab(iat)%gylm
675 !  Pawfgrtab(iat)%gylmgr
676 !  Pawfgrtab(iat)%gylmgr2
677 !  Pawfgrtab(ia)%nhatfr
678 !  Pawfgrtab(ia)%nhatfrgr
679 !  Pawfgrtab(ia)%expiqr
680  end do
681 
682 !Destroy atom table
683  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
684 
685 end subroutine pawfgrtab_print

m_pawfgrtab/pawfgrtab_redistribute [ Functions ]

[ Top ] [ m_pawfgrtab ] [ Functions ]

NAME

 pawfgrtab_redistribute

FUNCTION

   Redistribute an array of pawfgrtab datastructures
   Input pawfgrtab is given on a MPI communicator
   Output pawfgrtab is redistributed on another MPI communicator

INPUTS

  mpi_comm_in= input MPI (atom) communicator
  mpi_comm_out= output MPI (atom) communicator
  mpi_atmtab_in= --optional-- indexes of the input pawfgrtab treated by current proc
                 if not present, will be calculated in the present routine
  mpi_atmtab_out= --optional-- indexes of the output pawfgrtab treated by current proc
                  if not present, will be calculated in the present routine
  natom= --optional-- total number of atoms
  ----- Optional arguments used only for asynchronous communications -----
    RecvAtomProc(:)= rank of processor from which I expect atom (in mpi_comm_in)
    RecvAtomList(:)= indexes of atoms to be received by me
      RecvAtomList(irecv) are the atoms I expect from RecvAtomProc(irecv)
    SendAtomProc(:)= ranks of process destination of atom (in mpi_comm_in)
    SendAtomList(:)= indexes of atoms to be sent by me
      SendAtomList(isend) are the atoms sent to SendAtomProc(isend)

OUTPUT

  [pawfgrtab_out(:)]<type(pawfgrtab_type)>= --optional--
                    if present, the redistributed datastructure does not replace
                    the input one but is delivered in pawfgrtab_out
                    if not present, input and output datastructure are the same.

SIDE EFFECTS

  pawfgrtab(:)<type(pawfgrtab_type)>= input (and eventually output) pawfgrtab datastructures

PARENTS

      m_paral_pert

CHILDREN

SOURCE

1051 subroutine pawfgrtab_redistribute(pawfgrtab,mpi_comm_in,mpi_comm_out,&
1052 &                    natom,mpi_atmtab_in,mpi_atmtab_out,pawfgrtab_out,&
1053 &                    SendAtomProc,SendAtomList,RecvAtomProc,RecvAtomList)
1054 
1055 
1056 !This section has been created automatically by the script Abilint (TD).
1057 !Do not modify the following lines by hand.
1058 #undef ABI_FUNC
1059 #define ABI_FUNC 'pawfgrtab_redistribute'
1060 !End of the abilint section
1061 
1062  implicit none
1063 
1064 !Arguments ------------------------------------
1065 !scalars
1066  integer,intent(in) :: mpi_comm_in,mpi_comm_out
1067  integer,optional,intent(in) :: natom
1068 !arrays
1069  integer,intent(in),optional,target :: mpi_atmtab_in(:),mpi_atmtab_out(:)
1070  integer,intent(in),optional :: SendAtomProc(:),SendAtomList(:),RecvAtomProc(:),RecvAtomList(:)
1071  type(pawfgrtab_type),allocatable,intent(inout) :: pawfgrtab(:)
1072  type(pawfgrtab_type),pointer,optional :: pawfgrtab_out(:)
1073 
1074 !Local variables-------------------------------
1075 !scalars
1076  integer :: algo_option,iat_in,iat_out,iatom,i1,ierr,iircv,iisend,imsg,imsg1,imsg_current
1077  integer :: iproc_rcv,iproc_send,ireq,me_exch,mpi_comm_exch,my_natom_in,my_natom_out,my_tag,natom_tot,nb_dp
1078  integer :: nb_int,nb_msg,nbmsg_incoming,nbrecvmsg,nbsendreq,nbsend,nbsent,nbrecv,next,npawfgrtab_sent
1079  integer :: nproc_in,nproc_out
1080  logical :: flag,in_place,message_yet_prepared,my_atmtab_in_allocated,my_atmtab_out_allocated,paral_atom
1081 !arrays
1082  integer :: buf_size(3),request1(3)
1083  integer,pointer :: my_atmtab_in(:),my_atmtab_out(:)
1084  integer,allocatable :: atmtab_send(:),atm_indx_in(:),atm_indx_out(:),buf_int1(:),From(:),request(:)
1085  integer,allocatable,target:: buf_int(:)
1086  integer,pointer :: buf_ints(:)
1087  logical,allocatable :: msg_pick(:)
1088  real(dp),allocatable :: buf_dp1(:)
1089  real(dp),allocatable,target :: buf_dp(:)
1090  real(dp),pointer :: buf_dps(:)
1091  type(coeffi1_type),target,allocatable :: tab_buf_int(:),tab_buf_atom(:)
1092  type(coeff1_type),target,allocatable :: tab_buf_dp(:)
1093  type(pawfgrtab_type),allocatable :: pawfgrtab_all(:)
1094  type(pawfgrtab_type),pointer :: pawfgrtab_out1(:)
1095 
1096 ! *************************************************************************
1097 
1098 !@pawfgrtab_type
1099 
1100  in_place=(.not.present(pawfgrtab_out))
1101  my_natom_in=size(pawfgrtab)
1102 
1103 !If not "in_place", destroy ) then
1104  if (.not.in_place) then
1105    if (associated(pawfgrtab_out)) then
1106      call pawfgrtab_free(pawfgrtab_out)
1107      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_out)
1108    end if
1109  end if
1110 
1111 !Special sequential case
1112  if (mpi_comm_in==xmpi_comm_self.and.mpi_comm_out==xmpi_comm_self) then
1113    if ((.not.in_place).and.(my_natom_in>0)) then
1114      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_in))
1115      call pawfgrtab_nullify(pawfgrtab_out)
1116      call pawfgrtab_copy(pawfgrtab,pawfgrtab_out)
1117    end if
1118    return
1119  end if
1120 
1121 !Get total natom
1122  if (present(natom)) then
1123    natom_tot=natom
1124  else
1125    natom_tot=my_natom_in
1126    call xmpi_sum(natom_tot,mpi_comm_in,ierr)
1127  end if
1128 
1129 !Select input distribution
1130  if (present(mpi_atmtab_in)) then
1131    my_atmtab_in => mpi_atmtab_in
1132    my_atmtab_in_allocated=.false.
1133  else
1134    call get_my_atmtab(mpi_comm_in,my_atmtab_in,my_atmtab_in_allocated,&
1135 &                     paral_atom,natom_tot,my_natom_in)
1136  end if
1137 
1138 !Select output distribution
1139  if (present(mpi_atmtab_out)) then
1140    my_natom_out=size(mpi_atmtab_out)
1141    my_atmtab_out => mpi_atmtab_out
1142    my_atmtab_out_allocated=.false.
1143  else
1144    call get_my_natom(mpi_comm_out,my_natom_out,natom_tot)
1145    call get_my_atmtab(mpi_comm_out,my_atmtab_out,my_atmtab_out_allocated,&
1146 &                     paral_atom,natom_tot)
1147  end if
1148 
1149 !Select algo according to optional input arguments
1150  algo_option=1
1151  if (present(SendAtomProc).and.present(SendAtomList).and.&
1152 &    present(RecvAtomProc).and.present(RecvAtomList)) algo_option=2
1153 
1154 
1155 !Brute force algorithm (allgather + scatter)
1156 !---------------------------------------------------------
1157  if (algo_option==1) then
1158 
1159    LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_all,(natom_tot))
1160    call pawfgrtab_nullify(pawfgrtab_all)
1161    call pawfgrtab_copy(pawfgrtab,pawfgrtab_all,comm_atom=mpi_comm_in,mpi_atmtab=my_atmtab_in)
1162    if (in_place) then
1163     call pawfgrtab_free(pawfgrtab)
1164     LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab)
1165     LIBPAW_DATATYPE_ALLOCATE(pawfgrtab,(my_natom_out))
1166     call pawfgrtab_nullify(pawfgrtab)
1167     call pawfgrtab_copy(pawfgrtab_all,pawfgrtab,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out)
1168    else
1169      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_out))
1170      call pawfgrtab_nullify(pawfgrtab_out)
1171      call pawfgrtab_copy(pawfgrtab_all,pawfgrtab_out,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out)
1172    end if
1173    call pawfgrtab_free(pawfgrtab_all)
1174    LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_all)
1175 
1176 
1177 !Asynchronous algorithm (asynchronous communications)
1178 !---------------------------------------------------------
1179  else if (algo_option==2) then
1180 
1181    nbsend=size(SendAtomProc) ; nbrecv=size(RecvAtomProc)
1182 
1183    if (in_place) then
1184      if (my_natom_out > 0) then
1185        LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out1,(my_natom_out))
1186        call pawfgrtab_nullify(pawfgrtab_out1)
1187      else
1188        LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out1,(0))
1189      end if
1190    else
1191      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_out))
1192      call pawfgrtab_nullify(pawfgrtab_out)
1193      pawfgrtab_out1=>pawfgrtab_out
1194    end if
1195 
1196    nproc_in=xmpi_comm_size(mpi_comm_in)
1197    nproc_out=xmpi_comm_size(mpi_comm_out)
1198    if (nproc_in<=nproc_out) mpi_comm_exch=mpi_comm_out
1199    if (nproc_in>nproc_out) mpi_comm_exch=mpi_comm_in
1200    me_exch=xmpi_comm_rank(mpi_comm_exch)
1201 
1202 !  Dimension put to the maximum to send
1203    LIBPAW_ALLOCATE(atmtab_send,(nbsend))
1204    LIBPAW_ALLOCATE(atm_indx_in,(natom_tot))
1205    atm_indx_in=-1
1206    do iatom=1,my_natom_in
1207      atm_indx_in(my_atmtab_in(iatom))=iatom
1208    end do
1209    LIBPAW_ALLOCATE(atm_indx_out,(natom_tot))
1210    atm_indx_out=-1
1211    do iatom=1,my_natom_out
1212      atm_indx_out(my_atmtab_out(iatom))=iatom
1213    end do
1214 
1215    LIBPAW_DATATYPE_ALLOCATE(tab_buf_int,(nbsend))
1216    LIBPAW_DATATYPE_ALLOCATE(tab_buf_dp,(nbsend))
1217    LIBPAW_DATATYPE_ALLOCATE(tab_buf_atom,(nbsend))
1218    LIBPAW_ALLOCATE(request,(3*nbsend))
1219 
1220 !  A send buffer in an asynchrone communication couldn't be deallocate before it has been receive
1221    nbsent=0 ; ireq=0 ; iisend=0 ; nbsendreq=0 ; nb_msg=0
1222    do iisend=1,nbsend
1223      iproc_rcv=SendAtomProc(iisend)
1224      next=-1
1225      if (iisend < nbsend) next=SendAtomProc(iisend+1)
1226      if (iproc_rcv /= me_exch) then
1227        nbsent=nbsent+1
1228        atmtab_send(nbsent)=SendAtomList(iisend) ! we groups the atoms sends to the same process
1229        if (iproc_rcv /= next) then
1230          if (nbsent > 0) then
1231 !          Check if message has been yet prepared
1232            message_yet_prepared=.false.
1233            do imsg=1,nb_msg
1234              if (size(tab_buf_atom(imsg)%value) /= nbsent) then
1235                cycle
1236              else
1237                do imsg1=1,nbsent
1238                  if (tab_buf_atom(imsg)%value(imsg1)/=atmtab_send(imsg1)) exit
1239                  message_yet_prepared=.true.
1240                  imsg_current=imsg
1241                end do
1242              end if
1243            enddo
1244 !          Create the message
1245            if (.not.message_yet_prepared) then
1246              nb_msg=nb_msg+1
1247              call pawfgrtab_isendreceive_fillbuffer( &
1248 &                   pawfgrtab,atmtab_send,atm_indx_in,nbsent,buf_int,nb_int,buf_dp,nb_dp)
1249              LIBPAW_ALLOCATE(tab_buf_int(nb_msg)%value,(nb_int))
1250              LIBPAW_ALLOCATE(tab_buf_dp(nb_msg)%value,(nb_dp))
1251              tab_buf_int(nb_msg)%value(1:nb_int)=buf_int(1:nb_int)
1252              tab_buf_dp(nb_msg)%value(1:nb_dp)=buf_dp(1:nb_dp)
1253              LIBPAW_DEALLOCATE(buf_int)
1254              LIBPAW_DEALLOCATE(buf_dp)
1255              LIBPAW_ALLOCATE(tab_buf_atom(nb_msg)%value, (nbsent))
1256              tab_buf_atom(nb_msg)%value(1:nbsent)=atmtab_send(1:nbsent)
1257              imsg_current=nb_msg
1258            end if
1259 !          Communicate
1260            buf_size(1)=size(tab_buf_int(imsg_current)%value)
1261            buf_size(2)=size(tab_buf_dp(imsg_current)%value)
1262            buf_size(3)=nbsent
1263            buf_ints=>tab_buf_int(imsg_current)%value
1264            buf_dps=>tab_buf_dp(imsg_current)%value
1265            my_tag=400
1266            ireq=ireq+1
1267            call xmpi_isend(buf_size,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1268            my_tag=401
1269            ireq=ireq+1
1270            call xmpi_isend(buf_ints,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1271            my_tag=402
1272            ireq=ireq+1
1273            call xmpi_isend(buf_dps,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1274            nbsendreq=ireq
1275            nbsent=0
1276          end if
1277        end if
1278      else ! Just a renumbering, not a sending
1279        iat_in=atm_indx_in(SendAtomList(iisend))
1280        iat_out=atm_indx_out(my_atmtab_in(iat_in))
1281        call pawfgrtab_copy(pawfgrtab(iat_in:iat_in),pawfgrtab_out1(iat_out:iat_out))
1282        nbsent=0
1283      end if
1284    end do
1285 
1286    LIBPAW_ALLOCATE(From,(nbrecv))
1287    From(:)=-1 ; nbrecvmsg=0
1288    do iircv=1,nbrecv
1289      iproc_send=RecvAtomProc(iircv) !receive from (RcvAtomProc is sorted by growing process)
1290      next=-1
1291      if (iircv < nbrecv) next=RecvAtomProc(iircv+1)
1292      if (iproc_send /= me_exch .and. iproc_send/=next) then
1293        nbrecvmsg=nbrecvmsg+1
1294        From(nbrecvmsg)=iproc_send
1295      end if
1296    end do
1297 
1298    LIBPAW_ALLOCATE(msg_pick,(nbrecvmsg))
1299    msg_pick=.false.
1300    nbmsg_incoming=nbrecvmsg
1301    do while (nbmsg_incoming > 0)
1302      do i1=1,nbrecvmsg
1303        if (.not.msg_pick(i1)) then
1304          iproc_send=From(i1)
1305          flag=.false.
1306          my_tag=400
1307          call xmpi_iprobe(iproc_send,my_tag,mpi_comm_exch,flag,ierr)
1308          if (flag) then
1309            msg_pick(i1)=.true.
1310            call xmpi_irecv(buf_size,iproc_send,my_tag,mpi_comm_exch,request1(1),ierr)
1311            call xmpi_wait(request1(1),ierr)
1312            nb_int=buf_size(1)
1313            nb_dp=buf_size(2)
1314            npawfgrtab_sent=buf_size(3)
1315            LIBPAW_ALLOCATE(buf_int1,(nb_int))
1316            LIBPAW_ALLOCATE(buf_dp1,(nb_dp))
1317            my_tag=401
1318            call xmpi_irecv(buf_int1,iproc_send,my_tag,mpi_comm_exch,request1(2),ierr)
1319            my_tag=402
1320            call xmpi_irecv(buf_dp1,iproc_send,my_tag,mpi_comm_exch,request1(3),ierr)
1321            call xmpi_waitall(request1(2:3),ierr)
1322            call pawfgrtab_isendreceive_getbuffer(pawfgrtab_out1,npawfgrtab_sent,atm_indx_out, buf_int1,buf_dp1)
1323            nbmsg_incoming=nbmsg_incoming-1
1324            LIBPAW_DEALLOCATE(buf_int1)
1325            LIBPAW_DEALLOCATE(buf_dp1)
1326          end if
1327        end if
1328      end do
1329    end do
1330    LIBPAW_DEALLOCATE(msg_pick)
1331 
1332    if (in_place) then
1333      call pawfgrtab_free(pawfgrtab)
1334      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab)
1335      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab,(my_natom_out))
1336      call pawfgrtab_nullify(pawfgrtab)
1337      call pawfgrtab_copy(pawfgrtab_out1,pawfgrtab)
1338      call pawfgrtab_free(pawfgrtab_out1)
1339      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_out1)
1340    end if
1341 
1342 !  Wait for deallocating arrays that all sending operations has been realized
1343    if (nbsendreq > 0) then
1344      call xmpi_waitall(request(1:nbsendreq),ierr)
1345    end if
1346 
1347 !  Deallocate buffers
1348    do i1=1,nb_msg
1349      LIBPAW_DEALLOCATE(tab_buf_int(i1)%value)
1350      LIBPAW_DEALLOCATE(tab_buf_dp(i1)%value)
1351      LIBPAW_DEALLOCATE(tab_buf_atom(i1)%value)
1352    end do
1353    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_int)
1354    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_dp)
1355    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_atom)
1356    LIBPAW_DEALLOCATE(From)
1357    LIBPAW_DEALLOCATE(request)
1358    LIBPAW_DEALLOCATE(atmtab_send)
1359    LIBPAW_DEALLOCATE(atm_indx_in)
1360    LIBPAW_DEALLOCATE(atm_indx_out)
1361 
1362  end if !algo_option
1363 
1364 !Eventually release temporary pointers
1365  call free_my_atmtab(my_atmtab_in,my_atmtab_in_allocated)
1366  call free_my_atmtab(my_atmtab_out,my_atmtab_out_allocated)
1367 
1368 end subroutine pawfgrtab_redistribute

m_pawfgrtab/pawfgrtab_type [ Types ]

[ Top ] [ m_pawfgrtab ] [ Types ]

NAME

 pawfgrtab_type

FUNCTION

 For PAW, various arrays giving data related to fine grid for a given atom

SOURCE

 52  type,public :: pawfgrtab_type
 53 
 54 !Integer scalars
 55 
 56   integer :: cplex
 57    ! cplex=1 if potentials/densities are real, 2 if they are complex
 58 
 59   integer :: expiqr_allocated
 60    ! 1 if expiqr() is allocated (and computed)
 61 
 62   integer :: itypat
 63    ! itypat=type of the atom
 64 
 65   integer :: l_size
 66    ! 1+maximum value of l leading to non zero Gaunt coeffs
 67    ! for the considered atom type
 68 
 69   integer :: gylm_allocated
 70    ! 1 if gylm() is allocated (and computed)
 71 
 72   integer :: gylmgr_allocated
 73    ! 1 if gylmgr() is allocated (and computed)
 74 
 75   integer :: gylmgr2_allocated
 76    ! 1 if gylmgr2() is allocated (and computed)
 77 
 78   integer :: nfgd
 79    ! Number of Fine rectangular GriD points
 80    ! in the paw sphere around considered atom
 81 
 82   integer :: nhatfr_allocated
 83    ! 1 if nhatfr() is allocated (and computed)
 84 
 85   integer :: nhatfrgr_allocated
 86    ! 1 if nhatfrgr() is allocated (and computed)
 87 
 88   integer :: nspden
 89    ! Number of spin-density components
 90 
 91   integer :: rfgd_allocated
 92    ! 1 if rfgd() is allocated (and computed)
 93 
 94 !Integer arrays
 95   !integer :: ngfftf(18)
 96   ! Info on the dense FFT mesh.
 97 
 98   integer, allocatable :: ifftsph(:)
 99    ! ifftsph(nfgd)
100    ! Array giving the FFT index (fine grid) of a point in the paw
101    ! sphere around considered atom (ifftsph=ix+n1*(iy-1+n2*(iz-1))
102 
103 !Real (real(dp)) arrays
104 
105   real(dp), allocatable :: expiqr(:,:)
106    ! expiqr(2,nfgd)
107    ! Gives exp(-i.q.r) on the fine rectangular grid
108    ! where q is the phonons wavevector
109 
110   real(dp), allocatable :: gylm(:,:)
111    ! gylm(nfgd,l_size*l_size)
112    ! Gives g_l(r)*Y_lm(r) on the fine rectangular grid
113    ! around considered atom
114 
115   real(dp), allocatable :: gylmgr(:,:,:)
116    ! gylmgr(3,nfgd,l_size*l_size)
117    ! Gives the gradient of g_l(r)*Y_lm(r) wrt cart. coordinates
118    ! on the fine rectangular grid around considered atom
119 
120   real(dp), allocatable :: gylmgr2(:,:,:)
121    ! gylmgr(6,nfgd,l_size*l_size)
122    ! Gives the second gradient of g_l(r)*Y_lm(r) wrt cart. coordinates
123    ! on the fine rectangular grid around considered atom
124 
125   real(dp), allocatable :: nhatfr(:,:)
126    ! nhatfr(nfgd,nspden)
127    ! Gives the "frozen part" of 1st-order compensation charge density
128    ! on the fine rectangular grid around considered atom
129    ! Only used in response function calculations
130 
131   real(dp), allocatable :: nhatfrgr(:,:,:)
132    ! nhatfrgr(3,nfgd,nspden)
133    ! Gives the gradients (wrt cart. coordinates)
134    ! of "frozen part" of 1st-order compensation charge density
135    ! on the fine rectangular grid around considered atom
136    ! Only used in response function calculations
137 
138   real(dp), allocatable :: rfgd(:,:)
139    ! r(3,nfgd)
140    ! Gives all R vectors (r-r_atom) on the Fine rectangular GriD
141    ! around considered atom in Cartesian coordinates.
142 
143  end type pawfgrtab_type
144 
145 !public procedures.
146  public :: pawfgrtab_init           ! Constructor
147  public :: pawfgrtab_free           ! Free memory
148  public :: pawfgrtab_nullify
149  public :: pawfgrtab_copy           ! Copy the object
150  public :: pawfgrtab_print          ! Print the content of the object.
151  public :: pawfgrtab_gather         ! MPI gather
152  public :: pawfgrtab_redistribute   ! MPI redistribution
153 
154 !private procedures.
155  private :: pawfgrtab_isendreceive_getbuffer
156  private :: pawfgrtab_isendreceive_fillbuffer