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-2022 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

23 #include "libpaw.h"
24 
25 MODULE m_pawfgrtab
26 
27  USE_DEFS
28  USE_MSG_HANDLING
29  USE_MPI_WRAPPERS
30  USE_MEMORY_PROFILING
31 
32  use m_paral_atom, only : get_my_atmtab, free_my_atmtab, get_my_natom
33 
34  implicit none
35 
36  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.

SOURCE

384 subroutine pawfgrtab_copy(pawfgrtab_in,pawfgrtab_cp, &
385 &                         mpi_atmtab,comm_atom) ! optional arguments
386 
387 !Arguments ------------------------------------
388 !scalars
389 integer,optional,intent(in) :: comm_atom
390 !arrays
391  integer,optional,target,intent(in) :: mpi_atmtab(:)
392  type(Pawfgrtab_type),intent(in) :: pawfgrtab_in(:)
393  type(Pawfgrtab_type),intent(inout),target :: pawfgrtab_cp(:)
394 
395 !Local variables-------------------------------
396 !scalars
397  integer :: ij,ij1,istat,l_size,l_size2,my_comm_atom,my_natom,nfgd,nspden,paral_case
398  integer :: siz_in, siz_max, siz_out
399  logical ::  my_atmtab_allocated,paral_atom
400  character(len=500) :: msg
401 !arrays
402  integer,pointer :: my_atmtab(:)
403  type(Pawfgrtab_type), pointer :: pawfgrtab_out(:)
404 
405 ! *************************************************************************
406 
407 !@Pawfgrtab_type
408 
409 !Retrieve sizes
410  siz_in=size(pawfgrtab_in);siz_out=size(pawfgrtab_cp)
411 
412 !Set up parallelism over atoms
413  paral_atom=(present(comm_atom));if (paral_atom) paral_atom=(xmpi_comm_size(comm_atom)>1)
414  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
415  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
416  my_atmtab_allocated=.false.
417 
418 !Determine in which case we are (parallelism, ...)
419 !No parallelism: a single copy operation
420  paral_case=0;siz_max=siz_in
421  pawfgrtab_out => pawfgrtab_cp
422  if (paral_atom) then
423    if (siz_out<siz_in) then ! Parallelism: the copy operation is a scatter
424      call get_my_natom(my_comm_atom,my_natom,siz_in)
425      if (my_natom==siz_out) then
426        call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,siz_in)
427        paral_case=1;siz_max=siz_out
428        pawfgrtab_out => pawfgrtab_cp
429      else
430        msg=' siz_out should be equal to my_natom !'
431        LIBPAW_BUG(msg)
432      end if
433    else                            ! Parallelism: the copy operation is a gather
434      call get_my_natom(my_comm_atom,my_natom,siz_out)
435      if (my_natom==siz_in) then
436        paral_case=2;siz_max=siz_in
437      else
438        msg=' siz_in should be equal to my_natom !'
439        LIBPAW_BUG(msg)
440      end if
441    end if
442  end if
443 
444 !First case: a simple copy or a scatter
445  if (siz_max > 0 .and. ((paral_case== 0).or.(paral_case==1))) then
446    call pawfgrtab_nullify(pawfgrtab_out)
447 
448 !  Loop on pawfgrtab components
449    do ij1=1,siz_max
450      ij=ij1; if (paral_case==1)ij=my_atmtab(ij1)
451 
452      pawfgrtab_out(ij1)%cplex=pawfgrtab_in(ij)%cplex
453      pawfgrtab_out(ij1)%expiqr_allocated=pawfgrtab_in(ij)%expiqr_allocated
454      pawfgrtab_out(ij1)%itypat=pawfgrtab_in(ij)%itypat
455      l_size=pawfgrtab_in(ij)%l_size
456      l_size2=l_size*l_size
457      pawfgrtab_out(ij1)%l_size=l_size
458      pawfgrtab_out(ij1)%gylm_allocated=pawfgrtab_in(ij)%gylm_allocated
459      pawfgrtab_out(ij1)%gylmgr_allocated=pawfgrtab_in(ij)%gylmgr_allocated
460      pawfgrtab_out(ij1)%gylmgr2_allocated=pawfgrtab_in(ij)%gylmgr2_allocated
461      pawfgrtab_out(ij1)%nfgd=pawfgrtab_in(ij)%nfgd
462      pawfgrtab_out(ij1)%nhatfr_allocated=pawfgrtab_in(ij)%nhatfr_allocated
463      pawfgrtab_out(ij1)%nhatfrgr_allocated=pawfgrtab_in(ij)%nhatfrgr_allocated
464      nspden=pawfgrtab_in(ij)%nspden
465      pawfgrtab_out(ij1)%nspden=nspden
466      pawfgrtab_out(ij1)%rfgd_allocated=pawfgrtab_in(ij)%rfgd_allocated
467      nfgd=pawfgrtab_in(ij)%nfgd
468      LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%ifftsph,(nfgd))
469      pawfgrtab_out(ij1)%ifftsph(:)=pawfgrtab_in(ij)%ifftsph(:)
470      if (pawfgrtab_out(ij1)%expiqr_allocated>0) then
471        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%expiqr,(2,nfgd))
472        pawfgrtab_out(ij1)%expiqr(:,:)=pawfgrtab_in(ij)%expiqr(:,:)
473      end if
474      if (pawfgrtab_out(ij1)%gylm_allocated>0) then
475        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylm,(nfgd,l_size2))
476        pawfgrtab_out(ij1)%gylm(:,:)=pawfgrtab_in(ij)%gylm(:,:)
477      end if
478      if (pawfgrtab_out(ij1)%gylmgr_allocated>0) then
479        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylmgr,(3,nfgd,l_size2))
480        pawfgrtab_out(ij1)%gylmgr(:,:,:)=pawfgrtab_in(ij)%gylmgr(:,:,:)
481      end if
482      if (pawfgrtab_out(ij1)%gylmgr2_allocated>0) then
483        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%gylmgr2,(6,nfgd,l_size2))
484        pawfgrtab_out(ij1)%gylmgr2(:,:,:)=pawfgrtab_in(ij)%gylmgr2(:,:,:)
485      end if
486      if (pawfgrtab_out(ij1)%nhatfrgr_allocated>0) then
487        LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%nhatfrgr,(3,nfgd,nspden))
488        pawfgrtab_out(ij1)%nhatfrgr(:,:,:)=pawfgrtab_in(ij)%nhatfrgr(:,:,:)
489       end if
490       if (pawfgrtab_out(ij1)%rfgd_allocated>0) then
491         LIBPAW_ALLOCATE(pawfgrtab_out(ij1)%rfgd,(3,nfgd))
492         pawfgrtab_out(ij1)%rfgd(:,:)=pawfgrtab_in(ij)%rfgd(:,:)
493       end if
494 
495     end do
496  end if
497 
498 !Second case: a gather
499  if (paral_case==2) then
500    call pawfgrtab_gather(pawfgrtab_in,pawfgrtab_cp,my_comm_atom,istat,my_atmtab)
501  end if
502 
503 !Destroy atom table used for parallelism
504  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
505 
506 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

SOURCE

261 subroutine pawfgrtab_free(Pawfgrtab)
262 
263 !Arguments ------------------------------------
264 !arrays
265  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
266 
267 !Local variables-------------------------------
268 !scalars
269  integer :: iat,natom
270 
271 ! *************************************************************************
272 
273 !@Pawfgrtab_type
274 
275  natom=SIZE(Pawfgrtab);if (natom==0) return
276  do iat=1,natom
277   if (allocated(Pawfgrtab(iat)%ifftsph))  then
278     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%ifftsph)
279   end if
280   if (allocated(Pawfgrtab(iat)%gylm   ))  then
281     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylm)
282   end if
283   if (allocated(Pawfgrtab(iat)%gylmgr ))  then
284     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylmgr)
285   end if
286   if (allocated(Pawfgrtab(iat)%gylmgr2))  then
287     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%gylmgr2)
288   end if
289   if (allocated(Pawfgrtab(iat)%nhatfr ))  then
290     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%nhatfr)
291   end if
292   if (allocated(Pawfgrtab(iat)%nhatfrgr))  then
293     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%nhatfrgr)
294   end if
295   if (allocated(Pawfgrtab(iat)%rfgd   ))  then
296     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%rfgd)
297   end if
298   if (allocated(Pawfgrtab(iat)%expiqr))   then
299     LIBPAW_DEALLOCATE(Pawfgrtab(iat)%expiqr)
300   end if
301   Pawfgrtab(iat)%gylm_allocated=0
302   Pawfgrtab(iat)%gylmgr_allocated=0
303   Pawfgrtab(iat)%gylmgr2_allocated=0
304   Pawfgrtab(iat)%nhatfr_allocated=0
305   Pawfgrtab(iat)%nhatfrgr_allocated=0
306   Pawfgrtab(iat)%rfgd_allocated=0
307   Pawfgrtab(iat)%expiqr_allocated=0
308  end do
309 
310 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

SOURCE

635 subroutine pawfgrtab_gather(pawfgrtab,pawfgrtab_gathered,comm_atom,istat, &
636 &                           mpi_atmtab) ! optional argument
637 
638 !Arguments ------------------------------------
639 !scalars
640  integer,intent(in) :: comm_atom
641  integer, intent(out) :: istat
642 !arrays
643  integer,optional,target,intent(in) :: mpi_atmtab(:)
644  type(pawfgrtab_type),intent(in) :: pawfgrtab(:)
645  type(pawfgrtab_type),intent(inout) :: pawfgrtab_gathered(:)
646 
647 !Local variables-------------------------------
648 !scalars
649  integer :: buf_int_size,buf_int_size_all,buf_dp_size,buf_dp_size_all,i1,i2,iat,iatot
650  integer :: ierr,ij,indx_int,indx_dp,l_size,l_size2,my_natom,natom,nfgd,nproc_atom,nspden
651  logical :: my_atmtab_allocated,paral_atom
652  character(len=500) :: msg
653 !arrays
654  integer :: bufsz(2)
655  integer,allocatable :: buf_int(:),buf_int_all(:)
656  integer,allocatable :: count_dp(:),count_int(:),count_tot(:),displ_dp(:),displ_int(:)
657  integer,pointer :: my_atmtab(:)
658  real(dp),allocatable :: buf_dp(:),buf_dp_all(:)
659 
660 ! *************************************************************************
661 
662 !@Pawfgrtab_type
663 
664  istat=0
665  my_natom=size(pawfgrtab);natom=size(pawfgrtab_gathered)
666 
667 !Set up parallelism over atoms
668  paral_atom=(my_natom/=natom)
669  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
670  my_atmtab_allocated = .False.
671  if (paral_atom) then
672    call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
673  end if
674  nproc_atom=xmpi_comm_size(comm_atom)
675 
676 !Without parallelism, just copy input to output
677  if (.not.paral_atom) then
678    istat=1;return
679  end if
680 
681 !Compute size of buffers
682  buf_int_size=0
683  do iat=1,my_natom
684    buf_int_size=buf_int_size+pawfgrtab(iat)%nfgd+13
685  end do
686  buf_dp_size=0
687  do iat=1,my_natom
688    if (pawfgrtab(iat)%expiqr_allocated>0) then
689      buf_dp_size=buf_dp_size + size(pawfgrtab(iat)%expiqr(:,:))
690    end if
691    if (pawfgrtab(iat)%gylm_allocated>0) then
692      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylm(:,:))
693    end if
694    if (pawfgrtab(iat)%gylmgr_allocated>0) then
695      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylmgr(:,:,:))
696    end if
697    if (pawfgrtab(iat)%gylmgr2_allocated>0) then
698      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%gylmgr2(:,:,:))
699    end if
700    if (pawfgrtab(iat)%nhatfr_allocated>0) then
701      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%nhatfr(:,:))
702    end if
703    if (pawfgrtab(iat)%nhatfrgr_allocated>0) then
704      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%nhatfrgr(:,:,:))
705    end if
706    if (pawfgrtab(iat)%rfgd_allocated>0) then
707      buf_dp_size=buf_dp_size+size(pawfgrtab(iat)%rfgd(:,:))
708    end if
709  end do
710 
711 !Fill in input buffers
712  LIBPAW_ALLOCATE(buf_int,(buf_int_size))
713  LIBPAW_ALLOCATE(buf_dp,(buf_dp_size))
714  indx_int=1;indx_dp=1
715  do iat=1,my_natom
716    l_size=pawfgrtab(iat)%l_size
717    nfgd=pawfgrtab(iat)%nfgd
718    nspden=pawfgrtab(iat)%nspden
719    buf_int(indx_int)=my_atmtab(iat); indx_int=indx_int+1;
720    buf_int(indx_int)=pawfgrtab(iat)%cplex; indx_int=indx_int+1;
721    buf_int(indx_int)=pawfgrtab(iat)%itypat; indx_int=indx_int+1;
722    buf_int(indx_int)=pawfgrtab(iat)%expiqr_allocated; indx_int=indx_int+1;
723    buf_int(indx_int)=pawfgrtab(iat)%l_size; indx_int=indx_int+1;
724    buf_int(indx_int)=pawfgrtab(iat)%gylm_allocated; indx_int=indx_int+1;
725    buf_int(indx_int)=pawfgrtab(iat)%gylmgr_allocated; indx_int=indx_int+1;
726    buf_int(indx_int)=pawfgrtab(iat)%gylmgr2_allocated; indx_int=indx_int+1;
727    buf_int(indx_int)=pawfgrtab(iat)%nfgd; indx_int=indx_int+1;
728    buf_int(indx_int)=pawfgrtab(iat)%nhatfr_allocated; indx_int=indx_int+1;
729    buf_int(indx_int)=pawfgrtab(iat)%nhatfrgr_allocated; indx_int=indx_int+1;
730    buf_int(indx_int)=pawfgrtab(iat)%nspden; indx_int=indx_int+1;
731    buf_int(indx_int)=pawfgrtab(iat)%rfgd_allocated; indx_int=indx_int+1;
732    if (nfgd>0) then
733      buf_int(indx_int:indx_int+nfgd-1)=pawfgrtab(iat)%ifftsph(1:nfgd)
734      indx_int=indx_int+nfgd;
735    end if
736    if (pawfgrtab(iat)%expiqr_allocated>0) then
737      do i1=1,nfgd
738        buf_dp(indx_dp:indx_dp+1)=pawfgrtab(iat)%expiqr(1:2,i1)
739        indx_dp=indx_dp+2
740      end do
741    end if
742    l_size2=l_size*l_size
743    if (pawfgrtab(iat)%gylm_allocated>0) then
744      do i1=1,l_size2
745        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab(iat)%gylm(1:nfgd,i1)
746       indx_dp=indx_dp+nfgd
747      end do
748    end if
749    if (pawfgrtab(iat)%gylmgr_allocated>0) then
750      do i1=1,l_size2
751        do i2=1,nfgd
752          buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%gylmgr(1:3,i2,i1)
753          indx_dp=indx_dp+3
754        end do
755      end do
756    end if
757    if (pawfgrtab(iat)%gylmgr2_allocated>0) then
758      do i1=1,l_size2
759        do i2=1,nfgd
760          buf_dp(indx_dp:indx_dp+5)=pawfgrtab(iat)%gylmgr2(1:6,i2,i1)
761          indx_dp=indx_dp+6
762        end do
763      end do
764    end if
765    if (pawfgrtab(iat)%nhatfr_allocated>0) then
766      do i1=1,nspden
767        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab(iat)%nhatfr(1:nfgd,i1)
768        indx_dp=indx_dp+nfgd
769      end do
770    end if
771    if (pawfgrtab(iat)%nhatfrgr_allocated>0) then
772      do i1=1,nspden
773        do i2=1,nfgd
774          buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%nhatfrgr(1:3,i2,i1)
775          indx_dp=indx_dp+3
776        end do
777      end do
778    end if
779    if (pawfgrtab(iat)%rfgd_allocated>0) then
780      do i1=1,nfgd
781        buf_dp(indx_dp:indx_dp+2)=pawfgrtab(iat)%rfgd(1:3,i1)
782        indx_dp=indx_dp+3
783      end do
784    end if
785  end do
786  if (indx_int/=1+buf_int_size) then
787    msg='Error (1) in pawfgrtab_gather: wrong buffer sizes !'
788    LIBPAW_BUG(msg)
789  end if
790  if (indx_dp/=1+buf_dp_size) then
791    msg='Error (2) in pawfgrtab_gather: wrong buffer sizes !'
792    LIBPAW_BUG(msg)
793  end if
794 
795 !Communicate (1 gather for integers, 1 gather for reals)
796  LIBPAW_ALLOCATE(count_int,(nproc_atom))
797  LIBPAW_ALLOCATE(displ_int,(nproc_atom))
798  LIBPAW_ALLOCATE(count_dp ,(nproc_atom))
799  LIBPAW_ALLOCATE(displ_dp ,(nproc_atom))
800  LIBPAW_ALLOCATE(count_tot,(2*nproc_atom))
801  bufsz(1)=buf_int_size; bufsz(2)=buf_dp_size
802  call xmpi_allgather(bufsz,2,count_tot,comm_atom,ierr)
803  do ij=1,nproc_atom
804    count_int(ij)=count_tot(2*ij-1)
805    count_dp (ij)=count_tot(2*ij)
806  end do
807  displ_int(1)=0;displ_dp(1)=0
808  do ij=2,nproc_atom
809    displ_int(ij)=displ_int(ij-1)+count_int(ij-1)
810    displ_dp (ij)=displ_dp (ij-1)+count_dp (ij-1)
811  end do
812  buf_int_size_all=sum(count_int)
813  buf_dp_size_all =sum(count_dp)
814  LIBPAW_DEALLOCATE(count_tot)
815  LIBPAW_ALLOCATE(buf_int_all,(buf_int_size_all))
816  LIBPAW_ALLOCATE(buf_dp_all ,(buf_dp_size_all))
817  call xmpi_allgatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,comm_atom,ierr)
818  call xmpi_allgatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,comm_atom,ierr)
819  LIBPAW_DEALLOCATE(count_int)
820  LIBPAW_DEALLOCATE(displ_int)
821  LIBPAW_DEALLOCATE(count_dp)
822  LIBPAW_DEALLOCATE(displ_dp)
823 
824 !Retrieve output datastructures
825  indx_int=1; indx_dp=1
826  call pawfgrtab_free(pawfgrtab_gathered)
827  call pawfgrtab_nullify(pawfgrtab_gathered)
828  do iat=1,natom
829    iatot=buf_int_all(indx_int); indx_int=indx_int+1
830    pawfgrtab_gathered(iatot)%cplex=buf_int_all(indx_int); indx_int=indx_int+1;
831    pawfgrtab_gathered(iatot)%itypat=buf_int_all(indx_int); indx_int=indx_int+1;
832    pawfgrtab_gathered(iatot)%expiqr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
833    pawfgrtab_gathered(iatot)%l_size=buf_int_all(indx_int); indx_int=indx_int+1;
834    pawfgrtab_gathered(iatot)%gylm_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
835    pawfgrtab_gathered(iatot)%gylmgr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
836    pawfgrtab_gathered(iatot)%gylmgr2_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
837    pawfgrtab_gathered(iatot)%nfgd=buf_int_all(indx_int); indx_int=indx_int+1;
838    pawfgrtab_gathered(iatot)%nhatfr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
839    pawfgrtab_gathered(iatot)%nhatfrgr_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
840    pawfgrtab_gathered(iatot)%nspden=buf_int_all(indx_int); indx_int=indx_int+1;
841    pawfgrtab_gathered(iatot)%rfgd_allocated=buf_int_all(indx_int); indx_int=indx_int+1;
842    nfgd=pawfgrtab_gathered(iatot)%nfgd
843    l_size=pawfgrtab_gathered(iatot)%l_size
844    nspden= pawfgrtab_gathered(iatot)%nspden
845    l_size2=l_size*l_size
846    if (nfgd>0) then
847      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%ifftsph,(nfgd))
848      pawfgrtab_gathered(iatot)%ifftsph(1:nfgd)=buf_int_all(indx_int:indx_int+nfgd-1)
849      indx_int=indx_int+nfgd
850    end if
851    if (pawfgrtab_gathered(iatot)%expiqr_allocated>0) then
852      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%expiqr,(2,nfgd))
853      do i1=1,nfgd
854        pawfgrtab_gathered(iatot)%expiqr(1:2,i1)= buf_dp_all(indx_dp:indx_dp+1)
855        indx_dp=indx_dp+2
856      end do
857    end if
858    if (pawfgrtab_gathered(iatot)%gylm_allocated>0) then
859      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylm,(nfgd,l_size2))
860      do i1=1,l_size2
861        pawfgrtab_gathered(iatot)%gylm(1:nfgd,i1)=buf_dp_all(indx_dp:indx_dp+nfgd-1)
862        indx_dp=indx_dp+nfgd
863      end do
864    end if
865    if (pawfgrtab_gathered(iatot)%gylmgr_allocated>0) then
866      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylmgr,(3,nfgd,l_size2))
867      do i1=1,l_size2
868        do i2=1,nfgd
869          pawfgrtab_gathered(iatot)%gylmgr(1:3,i2,i1)=buf_dp_all(indx_dp:indx_dp+2)
870          indx_dp=indx_dp+3
871        end do
872      end do
873    end if
874    if (pawfgrtab_gathered(iatot)%gylmgr2_allocated>0) then
875      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%gylmgr2,(6,nfgd,l_size2))
876      do i1=1,l_size2
877        do i2=1,nfgd
878          pawfgrtab_gathered(iatot)%gylmgr2(1:6,i2,i1)=buf_dp_all(indx_dp:indx_dp+5)
879          indx_dp=indx_dp+6
880        end do
881      end do
882    end if
883    if (pawfgrtab_gathered(iatot)%nhatfr_allocated>0) then
884      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%nhatfr,(nfgd,nspden))
885      do i1=1,nspden
886        pawfgrtab_gathered(iatot)%nhatfr(1:nfgd,i1)=buf_dp_all(indx_dp:indx_dp+nfgd-1)
887        indx_dp=indx_dp+nfgd
888      end do
889    end if
890    if (pawfgrtab_gathered(iatot)%nhatfrgr_allocated>0) then
891      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%nhatfrgr,(3,nfgd,nspden))
892      do i1=1,nspden
893        do i2=1,nfgd
894          pawfgrtab_gathered(iatot)%nhatfrgr(1:3,i2,i1)=buf_dp_all(indx_dp:indx_dp+2)
895          indx_dp=indx_dp+3
896        end do
897      end do
898    end if
899    if (pawfgrtab_gathered(iatot)%rfgd_allocated>0) then
900      LIBPAW_ALLOCATE(pawfgrtab_gathered(iatot)%rfgd,(3,nfgd))
901      do i1=1,nfgd
902        pawfgrtab_gathered(iatot)%rfgd(1:3,i1)=buf_dp_all(indx_dp:indx_dp+2)
903        indx_dp=indx_dp+3
904      end do
905    end if
906  end do
907 
908 !Free memory
909  LIBPAW_DEALLOCATE(buf_int)
910  LIBPAW_DEALLOCATE(buf_int_all)
911  LIBPAW_DEALLOCATE(buf_dp)
912  LIBPAW_DEALLOCATE(buf_dp_all)
913 
914 !Destroy atom table
915  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
916 
917 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

SOURCE

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

SOURCE

1431 subroutine pawfgrtab_isendreceive_fillbuffer(pawfgrtab, atmtab_send,atm_indx_send,npawfgrtab_send,&
1432 &                                            buf_int,buf_int_size,buf_dp,buf_dp_size)
1433 
1434 !Arguments ------------------------------------
1435 !scalars
1436  integer,intent(out) :: buf_int_size,buf_dp_size
1437  integer,intent(in) :: npawfgrtab_send
1438 !arrays
1439  integer,intent(in) :: atmtab_send(:),atm_indx_send(:)
1440  integer,intent(out),allocatable :: buf_int(:)
1441  real(dp),intent(out),allocatable :: buf_dp(:)
1442  type(pawfgrtab_type),target,intent(inout) :: pawfgrtab(:)
1443 
1444 !Local variables-------------------------------
1445 !scalars
1446  integer :: i1,i2,iat,iat1,iatom_tot,indx_int,indx_dp,l_size,l_size2,nfgd,nspden
1447  character(len=500) :: msg
1448  type(pawfgrtab_type),pointer :: pawfgrtab1
1449 
1450 ! *********************************************************************
1451 
1452 !Compute size of buffers
1453  buf_int_size=0;buf_dp_size=0
1454  do iat1=1,npawfgrtab_send
1455    iatom_tot=atmtab_send(iat1)
1456    iat=atm_indx_send(iatom_tot)
1457    pawfgrtab1=>pawfgrtab(iat)
1458    buf_int_size=buf_int_size+pawfgrtab1%nfgd+13
1459    if (pawfgrtab1%expiqr_allocated>0) then
1460      buf_dp_size=buf_dp_size + size(pawfgrtab1%expiqr(:,:))
1461    end if
1462    if (pawfgrtab1%gylm_allocated>0) then
1463      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylm(:,:))
1464    end if
1465    if (pawfgrtab1%gylmgr_allocated>0) then
1466      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylmgr(:,:,:))
1467    end if
1468    if (pawfgrtab1%gylmgr2_allocated>0 ) then
1469      buf_dp_size=buf_dp_size+size(pawfgrtab1%gylmgr2(:,:,:))
1470    end if
1471    if (pawfgrtab1%nhatfr_allocated>0) then
1472      buf_dp_size=buf_dp_size+size(pawfgrtab1%nhatfr(:,:))
1473    end if
1474    if (pawfgrtab1%nhatfrgr_allocated>0) then
1475      buf_dp_size=buf_dp_size+size(pawfgrtab1%nhatfrgr(:,:,:))
1476    end if
1477    if (pawfgrtab1%rfgd_allocated>0) then
1478      buf_dp_size=buf_dp_size+size(pawfgrtab1%rfgd(:,:))
1479    end if
1480  end do
1481 
1482 !Fill in input buffers
1483  LIBPAW_ALLOCATE(buf_int,(buf_int_size))
1484  LIBPAW_ALLOCATE(buf_dp,(buf_dp_size))
1485  indx_int=1;indx_dp=1
1486  do iat1=1,npawfgrtab_send
1487     iatom_tot=atmtab_send(iat1)
1488    iat=atm_indx_send(iatom_tot)
1489    pawfgrtab1=>pawfgrtab(iat)
1490    l_size=pawfgrtab1%l_size
1491    nfgd=pawfgrtab1%nfgd
1492    nspden=pawfgrtab1%nspden
1493    buf_int(indx_int)=iatom_tot; indx_int=indx_int+1;
1494    buf_int(indx_int)=pawfgrtab1%cplex; indx_int=indx_int+1;
1495    buf_int(indx_int)=pawfgrtab1%itypat; indx_int=indx_int+1;
1496    buf_int(indx_int)=pawfgrtab1%expiqr_allocated; indx_int=indx_int+1;
1497    buf_int(indx_int)=pawfgrtab1%l_size; indx_int=indx_int+1;
1498    buf_int(indx_int)=pawfgrtab1%gylm_allocated; indx_int=indx_int+1;
1499    buf_int(indx_int)=pawfgrtab1%gylmgr_allocated; indx_int=indx_int+1;
1500    buf_int(indx_int)=pawfgrtab1%gylmgr2_allocated; indx_int=indx_int+1;
1501    buf_int(indx_int)=pawfgrtab1%nfgd; indx_int=indx_int+1;
1502    buf_int(indx_int)=pawfgrtab1%nhatfr_allocated; indx_int=indx_int+1;
1503    buf_int(indx_int)=pawfgrtab1%nhatfrgr_allocated; indx_int=indx_int+1;
1504    buf_int(indx_int)=pawfgrtab1%nspden; indx_int=indx_int+1;
1505    buf_int(indx_int)=pawfgrtab1%rfgd_allocated; indx_int=indx_int+1;
1506    if (nfgd>0) then
1507      buf_int(indx_int:indx_int+nfgd-1)=pawfgrtab1%ifftsph(1:nfgd)
1508      indx_int=indx_int+nfgd;
1509    end if
1510    if (pawfgrtab1%expiqr_allocated>0) then
1511      do i1=1,nfgd
1512        buf_dp(indx_dp:indx_dp+1)=pawfgrtab1%expiqr(1:2,i1)
1513        indx_dp=indx_dp+2
1514      end do
1515    end if
1516    l_size2=l_size*l_size
1517    if (pawfgrtab1%gylm_allocated>0) then
1518      do i1=1,l_size2
1519        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab1%gylm(1:nfgd,i1)
1520       indx_dp=indx_dp+nfgd
1521      end do
1522    end if
1523    if (pawfgrtab1%gylmgr_allocated>0) then
1524      do i1=1,l_size2
1525        do i2=1,nfgd
1526          buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%gylmgr(1:3,i2,i1)
1527          indx_dp=indx_dp+3
1528        end do
1529      end do
1530    end if
1531    if (pawfgrtab1%gylmgr2_allocated>0) then
1532      do i1=1,l_size2
1533        do i2=1,nfgd
1534          buf_dp(indx_dp:indx_dp+5)=pawfgrtab1%gylmgr2(1:6,i2,i1)
1535          indx_dp=indx_dp+6
1536        end do
1537      end do
1538    end if
1539    if (pawfgrtab1%nhatfr_allocated>0) then
1540      do i1=1,nspden
1541        buf_dp(indx_dp:indx_dp+nfgd-1)=pawfgrtab1%nhatfr(1:nfgd,i1)
1542        indx_dp=indx_dp+nfgd
1543      end do
1544    end if
1545    if (pawfgrtab1%nhatfrgr_allocated>0) then
1546      do i1=1,nspden
1547        do i2=1,nfgd
1548          buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%nhatfrgr(1:3,i2,i1)
1549          indx_dp=indx_dp+3
1550        end do
1551      end do
1552    end if
1553    if (pawfgrtab1%rfgd_allocated>0) then
1554      do i1=1,nfgd
1555        buf_dp(indx_dp:indx_dp+2)=pawfgrtab1%rfgd(1:3,i1)
1556        indx_dp=indx_dp+3
1557      end do
1558    end if
1559  end do
1560  if ((indx_int-1/=buf_int_size).or.(indx_dp-1/=buf_dp_size)) then
1561    write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size
1562    LIBPAW_BUG(msg)
1563  end if
1564 
1565 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

SOURCE

1293 subroutine pawfgrtab_isendreceive_getbuffer(pawfgrtab,npawfgrtab_send,atm_indx_recv,buf_int,buf_dp)
1294 
1295 !Arguments ------------------------------------
1296 !scalars
1297  integer,intent(in) :: npawfgrtab_send
1298 !arrays
1299  integer,intent(in) ::atm_indx_recv(:),buf_int(:)
1300  real(dp),intent(in) :: buf_dp(:)
1301  type(pawfgrtab_type),target,intent(inout) :: pawfgrtab(:)
1302 
1303 !Local variables-------------------------------
1304 !scalars
1305  integer :: buf_dp_size,buf_int_size,i1,i2,iat,iatom_tot,ij,indx_int,indx_dp
1306  integer :: l_size,l_size2,nfgd,nspden
1307  character(len=500) :: msg
1308  type(pawfgrtab_type),pointer :: pawfgrtab1
1309 
1310 ! *********************************************************************
1311 
1312  buf_int_size=size(buf_int)
1313  buf_dp_size=size(buf_dp)
1314  indx_int=1; indx_dp=1
1315 
1316  do ij=1,npawfgrtab_send
1317    iatom_tot=buf_int(indx_int); indx_int=indx_int+1
1318    iat= atm_indx_recv(iatom_tot)
1319    pawfgrtab1=>pawfgrtab(iat)
1320    pawfgrtab1%cplex=buf_int(indx_int); indx_int=indx_int+1;
1321    pawfgrtab1%itypat=buf_int(indx_int); indx_int=indx_int+1;
1322    pawfgrtab1%expiqr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1323    pawfgrtab1%l_size=buf_int(indx_int); indx_int=indx_int+1;
1324    pawfgrtab1%gylm_allocated=buf_int(indx_int); indx_int=indx_int+1;
1325    pawfgrtab1%gylmgr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1326    pawfgrtab1%gylmgr2_allocated=buf_int(indx_int); indx_int=indx_int+1;
1327    pawfgrtab1%nfgd=buf_int(indx_int); indx_int=indx_int+1;
1328    pawfgrtab1%nhatfr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1329    pawfgrtab1%nhatfrgr_allocated=buf_int(indx_int); indx_int=indx_int+1;
1330    pawfgrtab1%nspden=buf_int(indx_int); indx_int=indx_int+1;
1331    pawfgrtab1%rfgd_allocated=buf_int(indx_int); indx_int=indx_int+1;
1332    nfgd=pawfgrtab1%nfgd
1333    l_size=pawfgrtab1%l_size
1334    nspden= pawfgrtab1%nspden
1335    l_size2=l_size*l_size
1336    if (nfgd>0) then
1337      LIBPAW_ALLOCATE(pawfgrtab1%ifftsph,(nfgd))
1338      pawfgrtab1%ifftsph(1:nfgd)=buf_int(indx_int:indx_int+nfgd-1)
1339      indx_int=indx_int+nfgd
1340    end if
1341    if (pawfgrtab1%expiqr_allocated>0) then
1342      LIBPAW_ALLOCATE(pawfgrtab1%expiqr,(2,nfgd))
1343      do i1=1,nfgd
1344        pawfgrtab1%expiqr(1:2,i1)= buf_dp(indx_dp:indx_dp+1)
1345        indx_dp=indx_dp+2
1346      end do
1347    end if
1348    if (pawfgrtab1%gylm_allocated>0) then
1349      LIBPAW_ALLOCATE(pawfgrtab1%gylm,(nfgd,l_size2))
1350      do i1=1,l_size2
1351        pawfgrtab1%gylm(1:nfgd,i1)=buf_dp(indx_dp:indx_dp+nfgd-1)
1352        indx_dp=indx_dp+nfgd
1353      end do
1354    end if
1355    if (pawfgrtab1%gylmgr_allocated>0) then
1356      LIBPAW_ALLOCATE(pawfgrtab1%gylmgr,(3,nfgd,l_size2))
1357      do i1=1,l_size2
1358        do i2=1,nfgd
1359          pawfgrtab1%gylmgr(1:3,i2,i1)=buf_dp(indx_dp:indx_dp+2)
1360          indx_dp=indx_dp+3
1361        end do
1362      end do
1363    end if
1364    if (pawfgrtab1%gylmgr2_allocated>0) then
1365      LIBPAW_ALLOCATE(pawfgrtab1%gylmgr2,(6,nfgd,l_size2))
1366      do i1=1,l_size2
1367        do i2=1,nfgd
1368          pawfgrtab1%gylmgr2(1:6,i2,i1)=buf_dp(indx_dp:indx_dp+5)
1369          indx_dp=indx_dp+6
1370        end do
1371      end do
1372    end if
1373    if (pawfgrtab1%nhatfr_allocated>0) then
1374      LIBPAW_ALLOCATE(pawfgrtab1%nhatfr,(nfgd,nspden))
1375      do i1=1,nspden
1376        pawfgrtab1%nhatfr(1:nfgd,i1)=buf_dp(indx_dp:indx_dp+nfgd-1)
1377        indx_dp=indx_dp+nfgd
1378      end do
1379    end if
1380    if (pawfgrtab1%nhatfrgr_allocated>0) then
1381      LIBPAW_ALLOCATE(pawfgrtab1%nhatfrgr,(3,nfgd,nspden))
1382      do i1=1,nspden
1383        do i2=1,nfgd
1384          pawfgrtab1%nhatfrgr(1:3,i2,i1)=buf_dp(indx_dp:indx_dp+2)
1385          indx_dp=indx_dp+3
1386        end do
1387      end do
1388    end if
1389    if (pawfgrtab1%rfgd_allocated>0) then
1390      LIBPAW_ALLOCATE(pawfgrtab1%rfgd,(3,nfgd))
1391      do i1=1,nfgd
1392        pawfgrtab1%rfgd(1:3,i1)=buf_dp(indx_dp:indx_dp+2)
1393        indx_dp=indx_dp+3
1394      end do
1395    end if
1396  end do
1397  if ((indx_int/=1+buf_int_size).or.(indx_dp/=1+buf_dp_size)) then
1398    write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size
1399    LIBPAW_BUG(msg)
1400  end if
1401 
1402 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

SOURCE

327 subroutine pawfgrtab_nullify(Pawfgrtab)
328 
329 !Arguments ------------------------------------
330 !arrays
331  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
332 
333 !Local variables-------------------------------
334 !scalars
335  integer :: iat,natom
336 
337 ! *************************************************************************
338 
339 !@Pawfgrtab_type
340 
341 ! MGPAW: This one could be removed/renamed,
342 ! variables can be initialized in the datatype declaration
343 ! Do we need to expose this in the public API?
344 
345  natom=SIZE(Pawfgrtab);if (natom==0) return
346  do iat=1,natom
347   Pawfgrtab(iat)%nfgd              = 0
348   Pawfgrtab(iat)%gylm_allocated    = 0
349   Pawfgrtab(iat)%gylmgr_allocated  = 0
350   Pawfgrtab(iat)%gylmgr2_allocated = 0
351   Pawfgrtab(iat)%nhatfr_allocated  = 0
352   Pawfgrtab(iat)%nhatfrgr_allocated= 0
353   Pawfgrtab(iat)%rfgd_allocated    = 0
354   Pawfgrtab(iat)%expiqr_allocated  = 0
355  end do
356 
357 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)

SOURCE

533 subroutine pawfgrtab_print(Pawfgrtab,natom,unit,prtvol,mode_paral,mpi_atmtab,comm_atom)
534 
535 !Arguments ------------------------------------
536 !scalars
537  integer,optional,intent(in) :: comm_atom,natom,prtvol,unit
538  character(len=4),optional,intent(in) :: mode_paral
539 !arrays
540  type(Pawfgrtab_type),intent(inout) :: Pawfgrtab(:)
541  integer,optional,target,intent(in)  :: mpi_atmtab(:)
542 
543 !Local variables-------------------------------
544 !scalars
545  integer :: iat,iatom_tot,my_comm_atom,my_natom,my_unt,my_prtvol,size_Pawfgrtab
546  logical :: my_atmtab_allocated,paral_atom
547  character(len=4) :: my_mode
548  character(len=500) :: msg
549 !arrays
550  integer,pointer :: my_atmtab(:)
551 
552 ! *************************************************************************
553 
554 !@Pawfgrtab_type
555 
556  size_Pawfgrtab=SIZE(Pawfgrtab);if (size_Pawfgrtab==0) return
557 
558  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
559  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
560  my_mode  ='PERS' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
561  my_natom=size_Pawfgrtab; if (PRESENT(natom))  my_natom=natom
562 
563 !Set up parallelism over atoms
564  paral_atom=(present(comm_atom).and.my_natom/=size_Pawfgrtab)
565  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
566  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
567  call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,&
568 &                   my_natom,my_natom_ref=size_Pawfgrtab)
569 
570  write(msg,'(3a)')ch10,' === Content of the pawfgrtab datatype === ',ch10
571  call wrtout(my_unt,msg,my_mode)
572  do iat=1,my_natom
573    iatom_tot=iat;if(paral_atom) iatom_tot=my_atmtab(iat)
574 
575    write(msg,'(6(2a,i0))')ch10,&
576 &    ' > For atom number : ',iatom_tot,ch10,&
577 &    '    cplex= ',Pawfgrtab(iat)%cplex,ch10,&
578 &    '    # of spins components for density= ',Pawfgrtab(iat)%nspden,ch10,&
579 &    '    Type of atom= ',Pawfgrtab(iat)%itypat,ch10,&
580 &    '    1+ Max l in Gaunt coefficients= ',Pawfgrtab(iat)%l_size,ch10,&
581 &    '    Number of fine FFT points in PAW sphere= ',Pawfgrtab(iat)%nfgd
582    call wrtout(my_unt,msg,my_mode)
583 
584    if (my_prtvol>=3) then
585      write(msg,'(a,7(a,i2,a))')ch10,&
586 &      '    rfgd_allocated    : ',Pawfgrtab(iat)%rfgd_allocated,ch10,&
587 &      '    gylm_allocated    : ',Pawfgrtab(iat)%gylm_allocated,ch10,&
588 &      '    gylmgr_allocated  : ',Pawfgrtab(iat)%gylmgr_allocated,ch10,&
589 &      '    gylmgr2_allocated : ',Pawfgrtab(iat)%gylmgr2_allocated,ch10,&
590 &      '    nhatgr_allocated  : ',Pawfgrtab(iat)%nhatfr_allocated,ch10,&
591 &      '    nhatfrgr_allocated: ',Pawfgrtab(iat)%nhatfrgr_allocated,ch10,&
592 &      '    expiqr_allocated  : ',Pawfgrtab(iat)%expiqr_allocated,ch10
593      call wrtout(my_unt,msg,my_mode)
594 
595    end if
596 
597 !  These huge arrays are not printed out!
598 !  Pawfgrtab(iat)%ifftsph
599 !  Pawfgrtab(iat)%rfgd
600 !  Pawfgrtab(iat)%gylm
601 !  Pawfgrtab(iat)%gylmgr
602 !  Pawfgrtab(iat)%gylmgr2
603 !  Pawfgrtab(ia)%nhatfr
604 !  Pawfgrtab(ia)%nhatfrgr
605 !  Pawfgrtab(ia)%expiqr
606  end do
607 
608 !Destroy atom table
609  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
610 
611 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

SOURCE

 958 subroutine pawfgrtab_redistribute(pawfgrtab,mpi_comm_in,mpi_comm_out,&
 959 &                    natom,mpi_atmtab_in,mpi_atmtab_out,pawfgrtab_out,&
 960 &                    SendAtomProc,SendAtomList,RecvAtomProc,RecvAtomList)
 961 
 962 !Arguments ------------------------------------
 963 !scalars
 964  integer,intent(in) :: mpi_comm_in,mpi_comm_out
 965  integer,optional,intent(in) :: natom
 966 !arrays
 967  integer,intent(in),optional,target :: mpi_atmtab_in(:),mpi_atmtab_out(:)
 968  integer,intent(in),optional :: SendAtomProc(:),SendAtomList(:),RecvAtomProc(:),RecvAtomList(:)
 969  type(pawfgrtab_type),allocatable,intent(inout) :: pawfgrtab(:)
 970  type(pawfgrtab_type),pointer,optional :: pawfgrtab_out(:)
 971 
 972 !Local variables-------------------------------
 973 !scalars
 974  integer :: algo_option,iat_in,iat_out,iatom,i1,ierr,iircv,iisend,imsg,imsg1,imsg_current
 975  integer :: iproc_rcv,iproc_send,ireq,me_exch,mpi_comm_exch,my_natom_in,my_natom_out,my_tag,natom_tot,nb_dp
 976  integer :: nb_int,nb_msg,nbmsg_incoming,nbrecvmsg,nbsendreq,nbsend,nbsent,nbrecv,next,npawfgrtab_sent
 977  integer :: nproc_in,nproc_out
 978  logical :: flag,in_place,message_yet_prepared,my_atmtab_in_allocated,my_atmtab_out_allocated,paral_atom
 979 !arrays
 980  integer :: buf_size(3),request1(3)
 981  integer,pointer :: my_atmtab_in(:),my_atmtab_out(:)
 982  integer,allocatable :: atmtab_send(:),atm_indx_in(:),atm_indx_out(:),buf_int1(:),From(:),request(:)
 983  integer,allocatable,target:: buf_int(:)
 984  integer,pointer :: buf_ints(:)
 985  logical,allocatable :: msg_pick(:)
 986  real(dp),allocatable :: buf_dp1(:)
 987  real(dp),allocatable,target :: buf_dp(:)
 988  real(dp),pointer :: buf_dps(:)
 989  type(coeffi1_type),target,allocatable :: tab_buf_int(:),tab_buf_atom(:)
 990  type(coeff1_type),target,allocatable :: tab_buf_dp(:)
 991  type(pawfgrtab_type),allocatable :: pawfgrtab_all(:)
 992  type(pawfgrtab_type),pointer :: pawfgrtab_out1(:)
 993 
 994 ! *************************************************************************
 995 
 996 !@pawfgrtab_type
 997 
 998  in_place=(.not.present(pawfgrtab_out))
 999  my_natom_in=size(pawfgrtab)
1000 
1001 !If not "in_place", destroy ) then
1002  if (.not.in_place) then
1003    if (associated(pawfgrtab_out)) then
1004      call pawfgrtab_free(pawfgrtab_out)
1005      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_out)
1006    end if
1007  end if
1008 
1009 !Special sequential case
1010  if (mpi_comm_in==xmpi_comm_self.and.mpi_comm_out==xmpi_comm_self) then
1011    if ((.not.in_place).and.(my_natom_in>0)) then
1012      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_in))
1013      call pawfgrtab_nullify(pawfgrtab_out)
1014      call pawfgrtab_copy(pawfgrtab,pawfgrtab_out)
1015    end if
1016    return
1017  end if
1018 
1019 !Get total natom
1020  if (present(natom)) then
1021    natom_tot=natom
1022  else
1023    natom_tot=my_natom_in
1024    call xmpi_sum(natom_tot,mpi_comm_in,ierr)
1025  end if
1026 
1027 !Select input distribution
1028  if (present(mpi_atmtab_in)) then
1029    my_atmtab_in => mpi_atmtab_in
1030    my_atmtab_in_allocated=.false.
1031  else
1032    call get_my_atmtab(mpi_comm_in,my_atmtab_in,my_atmtab_in_allocated,&
1033 &                     paral_atom,natom_tot,my_natom_in)
1034  end if
1035 
1036 !Select output distribution
1037  if (present(mpi_atmtab_out)) then
1038    my_natom_out=size(mpi_atmtab_out)
1039    my_atmtab_out => mpi_atmtab_out
1040    my_atmtab_out_allocated=.false.
1041  else
1042    call get_my_natom(mpi_comm_out,my_natom_out,natom_tot)
1043    call get_my_atmtab(mpi_comm_out,my_atmtab_out,my_atmtab_out_allocated,&
1044 &                     paral_atom,natom_tot)
1045  end if
1046 
1047 !Select algo according to optional input arguments
1048  algo_option=1
1049  if (present(SendAtomProc).and.present(SendAtomList).and.&
1050 &    present(RecvAtomProc).and.present(RecvAtomList)) algo_option=2
1051 
1052 
1053 !Brute force algorithm (allgather + scatter)
1054 !---------------------------------------------------------
1055  if (algo_option==1) then
1056 
1057    LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_all,(natom_tot))
1058    call pawfgrtab_nullify(pawfgrtab_all)
1059    call pawfgrtab_copy(pawfgrtab,pawfgrtab_all,comm_atom=mpi_comm_in,mpi_atmtab=my_atmtab_in)
1060    if (in_place) then
1061     call pawfgrtab_free(pawfgrtab)
1062     LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab)
1063     LIBPAW_DATATYPE_ALLOCATE(pawfgrtab,(my_natom_out))
1064     call pawfgrtab_nullify(pawfgrtab)
1065     call pawfgrtab_copy(pawfgrtab_all,pawfgrtab,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out)
1066    else
1067      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_out))
1068      call pawfgrtab_nullify(pawfgrtab_out)
1069      call pawfgrtab_copy(pawfgrtab_all,pawfgrtab_out,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out)
1070    end if
1071    call pawfgrtab_free(pawfgrtab_all)
1072    LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_all)
1073 
1074 
1075 !Asynchronous algorithm (asynchronous communications)
1076 !---------------------------------------------------------
1077  else if (algo_option==2) then
1078 
1079    nbsend=size(SendAtomProc) ; nbrecv=size(RecvAtomProc)
1080 
1081    if (in_place) then
1082      if (my_natom_out > 0) then
1083        LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out1,(my_natom_out))
1084        call pawfgrtab_nullify(pawfgrtab_out1)
1085      else
1086        LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out1,(0))
1087      end if
1088    else
1089      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab_out,(my_natom_out))
1090      call pawfgrtab_nullify(pawfgrtab_out)
1091      pawfgrtab_out1=>pawfgrtab_out
1092    end if
1093 
1094    nproc_in=xmpi_comm_size(mpi_comm_in)
1095    nproc_out=xmpi_comm_size(mpi_comm_out)
1096    if (nproc_in<=nproc_out) mpi_comm_exch=mpi_comm_out
1097    if (nproc_in>nproc_out) mpi_comm_exch=mpi_comm_in
1098    me_exch=xmpi_comm_rank(mpi_comm_exch)
1099 
1100 !  Dimension put to the maximum to send
1101    LIBPAW_ALLOCATE(atmtab_send,(nbsend))
1102    LIBPAW_ALLOCATE(atm_indx_in,(natom_tot))
1103    atm_indx_in=-1
1104    do iatom=1,my_natom_in
1105      atm_indx_in(my_atmtab_in(iatom))=iatom
1106    end do
1107    LIBPAW_ALLOCATE(atm_indx_out,(natom_tot))
1108    atm_indx_out=-1
1109    do iatom=1,my_natom_out
1110      atm_indx_out(my_atmtab_out(iatom))=iatom
1111    end do
1112 
1113    LIBPAW_DATATYPE_ALLOCATE(tab_buf_int,(nbsend))
1114    LIBPAW_DATATYPE_ALLOCATE(tab_buf_dp,(nbsend))
1115    LIBPAW_DATATYPE_ALLOCATE(tab_buf_atom,(nbsend))
1116    LIBPAW_ALLOCATE(request,(3*nbsend))
1117 
1118 !  A send buffer in an asynchrone communication couldn't be deallocate before it has been receive
1119    nbsent=0 ; ireq=0 ; iisend=0 ; nbsendreq=0 ; nb_msg=0
1120    do iisend=1,nbsend
1121      iproc_rcv=SendAtomProc(iisend)
1122      next=-1
1123      if (iisend < nbsend) next=SendAtomProc(iisend+1)
1124      if (iproc_rcv /= me_exch) then
1125        nbsent=nbsent+1
1126        atmtab_send(nbsent)=SendAtomList(iisend) ! we groups the atoms sends to the same process
1127        if (iproc_rcv /= next) then
1128          if (nbsent > 0) then
1129 !          Check if message has been yet prepared
1130            message_yet_prepared=.false.
1131            do imsg=1,nb_msg
1132              if (size(tab_buf_atom(imsg)%value) /= nbsent) then
1133                cycle
1134              else
1135                do imsg1=1,nbsent
1136                  if (tab_buf_atom(imsg)%value(imsg1)/=atmtab_send(imsg1)) exit
1137                  message_yet_prepared=.true.
1138                  imsg_current=imsg
1139                end do
1140              end if
1141            enddo
1142 !          Create the message
1143            if (.not.message_yet_prepared) then
1144              nb_msg=nb_msg+1
1145              call pawfgrtab_isendreceive_fillbuffer( &
1146 &                   pawfgrtab,atmtab_send,atm_indx_in,nbsent,buf_int,nb_int,buf_dp,nb_dp)
1147              LIBPAW_ALLOCATE(tab_buf_int(nb_msg)%value,(nb_int))
1148              LIBPAW_ALLOCATE(tab_buf_dp(nb_msg)%value,(nb_dp))
1149              tab_buf_int(nb_msg)%value(1:nb_int)=buf_int(1:nb_int)
1150              tab_buf_dp(nb_msg)%value(1:nb_dp)=buf_dp(1:nb_dp)
1151              LIBPAW_DEALLOCATE(buf_int)
1152              LIBPAW_DEALLOCATE(buf_dp)
1153              LIBPAW_ALLOCATE(tab_buf_atom(nb_msg)%value, (nbsent))
1154              tab_buf_atom(nb_msg)%value(1:nbsent)=atmtab_send(1:nbsent)
1155              imsg_current=nb_msg
1156            end if
1157 !          Communicate
1158            buf_size(1)=size(tab_buf_int(imsg_current)%value)
1159            buf_size(2)=size(tab_buf_dp(imsg_current)%value)
1160            buf_size(3)=nbsent
1161            buf_ints=>tab_buf_int(imsg_current)%value
1162            buf_dps=>tab_buf_dp(imsg_current)%value
1163            my_tag=400
1164            ireq=ireq+1
1165            call xmpi_isend(buf_size,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1166            my_tag=401
1167            ireq=ireq+1
1168            call xmpi_isend(buf_ints,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1169            my_tag=402
1170            ireq=ireq+1
1171            call xmpi_isend(buf_dps,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr)
1172            nbsendreq=ireq
1173            nbsent=0
1174          end if
1175        end if
1176      else ! Just a renumbering, not a sending
1177        iat_in=atm_indx_in(SendAtomList(iisend))
1178        iat_out=atm_indx_out(my_atmtab_in(iat_in))
1179        call pawfgrtab_copy(pawfgrtab(iat_in:iat_in),pawfgrtab_out1(iat_out:iat_out))
1180        nbsent=0
1181      end if
1182    end do
1183 
1184    LIBPAW_ALLOCATE(From,(nbrecv))
1185    From(:)=-1 ; nbrecvmsg=0
1186    do iircv=1,nbrecv
1187      iproc_send=RecvAtomProc(iircv) !receive from (RcvAtomProc is sorted by growing process)
1188      next=-1
1189      if (iircv < nbrecv) next=RecvAtomProc(iircv+1)
1190      if (iproc_send /= me_exch .and. iproc_send/=next) then
1191        nbrecvmsg=nbrecvmsg+1
1192        From(nbrecvmsg)=iproc_send
1193      end if
1194    end do
1195 
1196    LIBPAW_ALLOCATE(msg_pick,(nbrecvmsg))
1197    msg_pick=.false.
1198    nbmsg_incoming=nbrecvmsg
1199    do while (nbmsg_incoming > 0)
1200      do i1=1,nbrecvmsg
1201        if (.not.msg_pick(i1)) then
1202          iproc_send=From(i1)
1203          flag=.false.
1204          my_tag=400
1205          call xmpi_iprobe(iproc_send,my_tag,mpi_comm_exch,flag,ierr)
1206          if (flag) then
1207            msg_pick(i1)=.true.
1208            call xmpi_irecv(buf_size,iproc_send,my_tag,mpi_comm_exch,request1(1),ierr)
1209            call xmpi_wait(request1(1),ierr)
1210            nb_int=buf_size(1)
1211            nb_dp=buf_size(2)
1212            npawfgrtab_sent=buf_size(3)
1213            LIBPAW_ALLOCATE(buf_int1,(nb_int))
1214            LIBPAW_ALLOCATE(buf_dp1,(nb_dp))
1215            my_tag=401
1216            call xmpi_irecv(buf_int1,iproc_send,my_tag,mpi_comm_exch,request1(2),ierr)
1217            my_tag=402
1218            call xmpi_irecv(buf_dp1,iproc_send,my_tag,mpi_comm_exch,request1(3),ierr)
1219            call xmpi_waitall(request1(2:3),ierr)
1220            call pawfgrtab_isendreceive_getbuffer(pawfgrtab_out1,npawfgrtab_sent,atm_indx_out, buf_int1,buf_dp1)
1221            nbmsg_incoming=nbmsg_incoming-1
1222            LIBPAW_DEALLOCATE(buf_int1)
1223            LIBPAW_DEALLOCATE(buf_dp1)
1224          end if
1225        end if
1226      end do
1227    end do
1228    LIBPAW_DEALLOCATE(msg_pick)
1229 
1230    if (in_place) then
1231      call pawfgrtab_free(pawfgrtab)
1232      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab)
1233      LIBPAW_DATATYPE_ALLOCATE(pawfgrtab,(my_natom_out))
1234      call pawfgrtab_nullify(pawfgrtab)
1235      call pawfgrtab_copy(pawfgrtab_out1,pawfgrtab)
1236      call pawfgrtab_free(pawfgrtab_out1)
1237      LIBPAW_DATATYPE_DEALLOCATE(pawfgrtab_out1)
1238    end if
1239 
1240 !  Wait for deallocating arrays that all sending operations has been realized
1241    if (nbsendreq > 0) then
1242      call xmpi_waitall(request(1:nbsendreq),ierr)
1243    end if
1244 
1245 !  Deallocate buffers
1246    do i1=1,nb_msg
1247      LIBPAW_DEALLOCATE(tab_buf_int(i1)%value)
1248      LIBPAW_DEALLOCATE(tab_buf_dp(i1)%value)
1249      LIBPAW_DEALLOCATE(tab_buf_atom(i1)%value)
1250    end do
1251    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_int)
1252    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_dp)
1253    LIBPAW_DATATYPE_DEALLOCATE(tab_buf_atom)
1254    LIBPAW_DEALLOCATE(From)
1255    LIBPAW_DEALLOCATE(request)
1256    LIBPAW_DEALLOCATE(atmtab_send)
1257    LIBPAW_DEALLOCATE(atm_indx_in)
1258    LIBPAW_DEALLOCATE(atm_indx_out)
1259 
1260  end if !algo_option
1261 
1262 !Eventually release temporary pointers
1263  call free_my_atmtab(my_atmtab_in,my_atmtab_in_allocated)
1264  call free_my_atmtab(my_atmtab_out,my_atmtab_out_allocated)
1265 
1266 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

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