TABLE OF CONTENTS


ABINIT/m_results_img [ Modules ]

[ Top ] [ Modules ]

NAME

  m_results_img

FUNCTION

  This module provides the definition of the results_img_type used
  to store results from GS calculations for a given image of the cell.

COPYRIGHT

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

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

26 #if defined HAVE_CONFIG_H
27 #include "config.h"
28 #endif
29 
30 #include "abi_common.h"
31 
32 MODULE m_results_img
33 
34  use defs_basis
35  use defs_abitypes
36  use m_energies
37  use m_abicore
38  use m_results_gs
39  use m_errors
40  use m_xmpi
41 
42  use m_geometry,   only : mkrdim, xred2xcart
43 
44 
45  implicit none
46 
47  private
48 
49 ! public procedures.
50  public :: init_results_img
51  public :: destroy_results_img
52  public :: nullify_results_img
53  public :: copy_results_img
54  public :: gather_results_img
55  public :: gather_array_img
56  public :: scatter_array_img
57  public :: get_geometry_img
58 
59  interface gather_array_img
60    module procedure gather_array_img_1D
61    module procedure gather_array_img_2D
62  end interface gather_array_img

m_results_img/copy_results_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  copy_results_img

FUNCTION

  Copy a results_img datastructure into another

INPUTS

  results_img_in=<type(results_img_type)>=input results_img datastructure

OUTPUT

  results_img_out=<type(results_img_type)>=output results_img datastructure

PARENTS

      gstateimg,m_results_img

CHILDREN

      mkrdim,xred2xcart

SOURCE

403 subroutine copy_results_img(results_img_in,results_img_out)
404 
405 
406 !This section has been created automatically by the script Abilint (TD).
407 !Do not modify the following lines by hand.
408 #undef ABI_FUNC
409 #define ABI_FUNC 'copy_results_img'
410 !End of the abilint section
411 
412  implicit none
413 
414 !Arguments ------------------------------------
415 !arrays
416  type(results_img_type),intent(in) :: results_img_in
417  type(results_img_type),intent(inout) :: results_img_out !vz_i
418 !Local variables-------------------------------
419 !scalars
420  integer :: natom_in,natom_out,npspalch_in,npspalch_out,ntypalch_in,ntypalch_out
421  integer :: nsppol_in, nsppol_out,ntypat_in,ntypat_out
422 
423 !************************************************************************
424 
425  !@results_img_type
426 
427  natom_in =results_img_in%natom
428  natom_out=results_img_out%natom
429  npspalch_in =results_img_in%npspalch
430  npspalch_out=results_img_out%npspalch
431  nsppol_in =results_img_in%nsppol
432  nsppol_out=results_img_out%nsppol
433  ntypalch_in =results_img_in%ntypalch
434  ntypalch_out=results_img_out%ntypalch
435  ntypat_in =results_img_in%ntypat
436  ntypat_out=results_img_out%ntypat
437 
438  if (natom_in>natom_out) then
439    if (allocated(results_img_out%xred))   then
440      ABI_DEALLOCATE(results_img_out%xred)
441    end if
442    if (allocated(results_img_out%vel))    then
443      ABI_DEALLOCATE(results_img_out%vel)
444    end if
445    if (allocated(results_img_out%vel_cell))    then
446      ABI_DEALLOCATE(results_img_out%vel_cell)
447    end if
448 
449    if (allocated(results_img_in%xred))   then
450      ABI_ALLOCATE(results_img_out%xred,(3,natom_in))
451    end if
452    if (allocated(results_img_in%vel))    then
453      ABI_ALLOCATE(results_img_out%vel,(3,natom_in))
454    end if
455    if (allocated(results_img_in%vel_cell))    then
456      ABI_ALLOCATE(results_img_out%vel_cell,(3,3))
457    end if
458  end if
459 
460  if (npspalch_in>npspalch_out .or. ntypalch_in>ntypalch_out) then
461    if (allocated(results_img_out%mixalch))   then
462      ABI_DEALLOCATE(results_img_out%mixalch)
463    end if
464  endif
465 
466  if (ntypat_in>ntypat_out) then
467    if (allocated(results_img_in%amu))    then
468      ABI_ALLOCATE(results_img_out%amu,(ntypat_in))
469    endif
470  endif
471 
472  results_img_out%natom  =results_img_in%natom
473  results_img_out%npspalch  =results_img_in%npspalch
474  results_img_out%nsppol =results_img_in%nsppol
475  results_img_out%ntypalch  =results_img_in%ntypalch
476  results_img_out%ntypat  =results_img_in%ntypat
477 
478  call copy_results_gs(results_img_in%results_gs,results_img_out%results_gs)
479 
480  if (allocated(results_img_in%acell)) results_img_out%acell(:)=results_img_in%acell(:)
481  if (allocated(results_img_in%amu))   results_img_out%amu(1:ntypat_in)=results_img_in%amu(1:ntypat_in)
482  if (allocated(results_img_in%mixalch)) &
483 &   results_img_out%mixalch(1:npspalch_in,1:ntypalch_in)=results_img_in%mixalch(1:npspalch_in,1:ntypalch_in)
484  if (allocated(results_img_in%rprim))   results_img_out%rprim(:,:)=results_img_in%rprim(:,:)
485  if (allocated(results_img_in%xred))    results_img_out%xred(:,1:natom_in)=results_img_in%xred(:,1:natom_in)
486  if (allocated(results_img_in%vel))     results_img_out%vel(:,1:natom_in)=results_img_in%vel(:,1:natom_in)
487  if (allocated(results_img_in%vel_cell))results_img_out%vel_cell(:,:)=results_img_in%vel_cell(:,:)
488 
489 end subroutine copy_results_img

m_results_img/destroy_results_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  destroy_results_img

FUNCTION

  Clean and destroy an array of results_img datastructures

INPUTS

OUTPUT

SIDE EFFECTS

  results_img(:)=<type(results_img_type)>=results_img datastructure array

PARENTS

      gstateimg,prtimg

CHILDREN

      mkrdim,xred2xcart

SOURCE

252 subroutine destroy_results_img(results_img)
253 
254 
255 !This section has been created automatically by the script Abilint (TD).
256 !Do not modify the following lines by hand.
257 #undef ABI_FUNC
258 #define ABI_FUNC 'destroy_results_img'
259 !End of the abilint section
260 
261  implicit none
262 
263 !Arguments ------------------------------------
264 !arrays
265  type(results_img_type),intent(inout) :: results_img(:)
266 !Local variables-------------------------------
267 !scalars
268  integer :: ii,results_img_size
269 
270 !************************************************************************
271 
272  !@results_img_type
273 
274  results_img_size=size(results_img)
275  if (results_img_size>0) then
276 
277    do ii=1,results_img_size
278 
279      results_img(ii)%natom=0
280      results_img(ii)%npspalch=0
281      results_img(ii)%nsppol  =0
282      results_img(ii)%ntypalch=0
283      results_img(ii)%ntypat=0
284 
285      if (associated(results_img(ii)%results_gs)) then
286        call destroy_results_gs(results_img(ii)%results_gs)
287        ABI_DATATYPE_DEALLOCATE(results_img(ii)%results_gs)
288      end if
289 
290      if (allocated(results_img(ii)%acell))  then
291        ABI_DEALLOCATE(results_img(ii)%acell)
292      end if
293      if (allocated(results_img(ii)%amu))  then
294        ABI_DEALLOCATE(results_img(ii)%amu)
295      end if
296      if (allocated(results_img(ii)%mixalch))  then
297        ABI_DEALLOCATE(results_img(ii)%mixalch)
298      end if
299      if (allocated(results_img(ii)%rprim))  then
300        ABI_DEALLOCATE(results_img(ii)%rprim)
301      end if
302      if (allocated(results_img(ii)%xred))   then
303        ABI_DEALLOCATE(results_img(ii)%xred)
304      end if
305      if (allocated(results_img(ii)%vel))    then
306        ABI_DEALLOCATE(results_img(ii)%vel)
307      end if
308      if (allocated(results_img(ii)%vel_cell))    then
309        ABI_DEALLOCATE(results_img(ii)%vel_cell)
310      end if
311    end do
312 
313  end if
314 
315 end subroutine destroy_results_img

m_results_img/gather_array_img_1D [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  gather_array_img_1D

FUNCTION

  Gather an real 1D-array (part of a results_img datastructure) using communicator
  over images (replicas) of the cell.
  Each contribution of single processor is gathered into a big array on master processor

INPUTS

  allgather= --optional, default=false--  if TRUE do ALL_GATHER instead of GATHER
  master= --optional, default=0-- index of master proc receiving gathered data (if allgather=false)
  mpi_enreg=information about MPI parallelization
  only_one_per_img= --optional, default=true--  if TRUE, the gather operation
                    is only done by one proc per image (master of the comm_cell)
  array_img(:,:)= (real) 1D-array distributed (has 2 dimensions; the 2nd one is nimage)

SIDE EFFECTS

  array_img_all(:,:)= (real) global (gathered) 1D-array
                      (has 2 dimensions; the 2nd one is nimagetot)

PARENTS

CHILDREN

      mkrdim,xred2xcart

SOURCE

812 subroutine gather_array_img_1D(array_img,array_img_all,mpi_enreg,&
813 &                              master,allgather,only_one_per_img) ! optional arguments
814 
815 
816 !This section has been created automatically by the script Abilint (TD).
817 !Do not modify the following lines by hand.
818 #undef ABI_FUNC
819 #define ABI_FUNC 'gather_array_img_1D'
820 !End of the abilint section
821 
822  implicit none
823 
824 !Arguments ------------------------------------
825 !scalars
826  integer,optional,intent(in) :: master
827  logical,optional,intent(in) :: allgather,only_one_per_img
828  type(MPI_type),intent(in) :: mpi_enreg
829 !arrays
830  real(dp),intent(in) :: array_img(:,:)
831  real(dp),intent(inout) :: array_img_all(:,:)
832 
833 !Local variables-------------------------------
834 !scalars
835  integer :: ibufr,ierr,iproc,jj
836  integer :: master_all,master_img,master_one_img,nimage,nimagetot
837  integer :: rsize,rsize_img,size1
838  logical :: do_allgather,i_am_master,one_per_img,use_array_all
839  !character(len=500) :: msg
840 !arrays
841  integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:)
842  real(dp),allocatable :: rbuffer(:),rbuffer_all(:)
843 
844 !************************************************************************
845 
846  !@results_img_type
847 
848  one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img
849  do_allgather=.false.;if (present(allgather)) do_allgather=allgather
850  master_all=0;if (present(master)) master_all=master
851  i_am_master=(mpi_enreg%me==master_all)
852 
853  master_img=0;master_one_img=0
854  use_array_all= &
855 &  (((     do_allgather).and.(     one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. &
856 &   ((     do_allgather).and.(.not.one_per_img))                                            .or. &
857 &   ((.not.do_allgather).and.(     one_per_img).and.(mpi_enreg%me==master_all))             .or. &
858 &   ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img)))
859 
860  size1=size(array_img,1)
861  if (use_array_all) then
862    if (size(array_img_all,1)/=size1) then
863      MSG_BUG('Wrong array_img_all size (1)')
864    end if
865  end if
866 
867  if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then
868 
869 !  Simple copy in case of 1 image
870    if (use_array_all) then
871      if (size(array_img_all,2)<=1) then
872        array_img_all(:,1)=array_img(:,1)
873        return
874      end if
875    endif
876 
877 !  Gather number of images treated by each proc
878    nimage=size(array_img,2)
879    ABI_ALLOCATE(nimage_all,(mpi_enreg%nproc_img))
880    call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr)
881    nimagetot=sum(nimage_all)
882    if (use_array_all) then
883      if (size(array_img_all,2)/=nimagetot) then
884        MSG_BUG('Wrong array_img_all size (2)!')
885      endif
886    end if
887 
888 !  Compute number of data
889    rsize=size1;rsize_img=nimage*rsize
890    ABI_ALLOCATE(rsize_img_all,(mpi_enreg%nproc_img))
891    rsize_img_all(:)=rsize*nimage_all(:)
892    ABI_DEALLOCATE(nimage_all)
893 
894 !  Compute shifts in buffer arrays for each proc
895    ABI_ALLOCATE(rbufshft,(mpi_enreg%nproc_img))
896    rbufshft(1)=0
897    do jj=2,mpi_enreg%nproc_img
898      rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1)
899    end do
900 
901 !  Load buffers
902    ABI_ALLOCATE(rbuffer,(rsize_img))
903    ibufr=0
904    do jj=1,nimage
905      rbuffer(ibufr+1:ibufr+rsize)=reshape(array_img(1:size1,jj),(/rsize/))
906      ibufr=ibufr+rsize
907    end do
908 
909 !  Gather all data
910    if (use_array_all)  then
911      ABI_ALLOCATE(rbuffer_all,(rsize*nimagetot))
912    end if
913    if (.not.use_array_all)  then
914      ABI_ALLOCATE(rbuffer_all,(0))
915    end if
916    if (do_allgather) then
917      call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
918 &                         mpi_enreg%comm_img,ierr)
919    else
920      call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
921 &                      master_img,mpi_enreg%comm_img,ierr)
922    end if
923    ABI_DEALLOCATE(rbuffer)
924    ABI_DEALLOCATE(rsize_img_all)
925 
926 !  Transfer buffers into gathered array_img_all (master proc only)
927    if (use_array_all) then
928      ABI_ALLOCATE(iimg,(mpi_enreg%nproc_img))
929      iimg=0
930      do jj=1,nimagetot
931 !      The following line supposes that images are sorted by increasing index
932        iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1
933        ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize
934        array_img_all(1:size1,jj)=reshape(rbuffer_all(ibufr+1:ibufr+rsize),(/size1/))
935      end do
936      ABI_DEALLOCATE(iimg)
937    end if
938 
939 !  Free memory
940    ABI_DEALLOCATE(rbufshft)
941    ABI_DEALLOCATE(rbuffer_all)
942 
943  end if
944 
945 end subroutine gather_array_img_1D

m_results_img/gather_array_img_2D [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  gather_array_img_2D

FUNCTION

  Gather an real 2D-array (part of a results_img datastructure) using communicator
  over images (replicas) of the cell.
  Each contribution of single processor is gathered into a big array on master processor

INPUTS

  allgather= --optional, default=false--  if TRUE do ALL_GATHER instead of GATHER
  master= --optional, default=0-- index of master proc receiving gathered data (if allgather=false)
  mpi_enreg=information about MPI parallelization
  only_one_per_img= --optional, default=true--  if TRUE, the gather operation
                    is only done by one proc per image (master of the comm_cell)
  array_img(:,:,:)= (real) 2D-array distributed (has 3 dimensions; the 3rd one is nimage)

SIDE EFFECTS

  array_img_all(:,:,:)= (real) global (gathered) 2D-array
                        (has 3 dimensions; the 3rd one is nimagetot)

PARENTS

CHILDREN

      mkrdim,xred2xcart

SOURCE

 978 subroutine gather_array_img_2D(array_img,array_img_all,mpi_enreg,&
 979 &                              master,allgather,only_one_per_img) ! optional arguments
 980 
 981 
 982 !This section has been created automatically by the script Abilint (TD).
 983 !Do not modify the following lines by hand.
 984 #undef ABI_FUNC
 985 #define ABI_FUNC 'gather_array_img_2D'
 986 !End of the abilint section
 987 
 988  implicit none
 989 
 990 !Arguments ------------------------------------
 991 !scalars
 992  integer,optional,intent(in) :: master
 993  logical,optional,intent(in) :: allgather,only_one_per_img
 994  type(MPI_type),intent(in) :: mpi_enreg
 995 !arrays
 996  real(dp),intent(in) :: array_img(:,:,:)
 997  real(dp),intent(inout) :: array_img_all(:,:,:)
 998 
 999 !Local variables-------------------------------
1000 !scalars
1001  integer :: ibufr,ierr,iproc,jj
1002  integer :: master_all,master_img,master_one_img,nimage,nimagetot
1003  integer :: rsize,rsize_img,size1,size2
1004  logical :: do_allgather,i_am_master,one_per_img,use_array_all
1005  !character(len=500) :: msg
1006 !arrays
1007  integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:)
1008  real(dp),allocatable :: rbuffer(:),rbuffer_all(:)
1009 
1010 !************************************************************************
1011 
1012  !@results_img_type
1013 
1014  one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img
1015  do_allgather=.false.;if (present(allgather)) do_allgather=allgather
1016  master_all=0;if (present(master)) master_all=master
1017  i_am_master=(mpi_enreg%me==master_all)
1018 
1019  master_img=0;master_one_img=0
1020  use_array_all= &
1021 &  (((     do_allgather).and.(     one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. &
1022 &   ((     do_allgather).and.(.not.one_per_img))                                            .or. &
1023 &   ((.not.do_allgather).and.(     one_per_img).and.(mpi_enreg%me==master_all))             .or. &
1024 &   ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img)))
1025 
1026  size1=size(array_img,1);size2=size(array_img,2)
1027  if (use_array_all) then
1028    if (size(array_img_all,1)/=size1.or.size(array_img_all,2)/=size2) then
1029      MSG_BUG('Wrong array_img_all size (1)!')
1030    end if
1031  end if
1032 
1033  if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then
1034 
1035 !  Simple copy in case of 1 image
1036    if (use_array_all) then
1037      if (size(array_img_all,3)<=1) then
1038        array_img_all(:,:,1)=array_img(:,:,1)
1039        return
1040      end if
1041    endif
1042 
1043 !  Gather number of images treated by each proc
1044    nimage=size(array_img,3)
1045    ABI_ALLOCATE(nimage_all,(mpi_enreg%nproc_img))
1046    call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr)
1047    nimagetot=sum(nimage_all)
1048    if (use_array_all) then
1049      if (size(array_img_all,3)/=nimagetot) then
1050        MSG_BUG('Wrong array_img_all size (2)!')
1051      endif
1052    end if
1053 
1054 !  Compute number of data
1055    rsize=size1*size2;rsize_img=nimage*rsize
1056    ABI_ALLOCATE(rsize_img_all,(mpi_enreg%nproc_img))
1057    rsize_img_all(:)=rsize*nimage_all(:)
1058    ABI_DEALLOCATE(nimage_all)
1059 
1060 !  Compute shifts in buffer arrays for each proc
1061    ABI_ALLOCATE(rbufshft,(mpi_enreg%nproc_img))
1062    rbufshft(1)=0
1063    do jj=2,mpi_enreg%nproc_img
1064      rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1)
1065    end do
1066 
1067 !  Load buffers
1068    ABI_ALLOCATE(rbuffer,(rsize_img))
1069    ibufr=0
1070    do jj=1,nimage
1071      rbuffer(ibufr+1:ibufr+rsize)=reshape(array_img(1:size1,1:size2,jj),(/rsize/))
1072      ibufr=ibufr+rsize
1073    end do
1074 
1075 !  Gather all data
1076    if (use_array_all)  then
1077      ABI_ALLOCATE(rbuffer_all,(rsize*nimagetot))
1078    end if
1079    if (.not.use_array_all)  then
1080      ABI_ALLOCATE(rbuffer_all,(0))
1081    end if
1082    if (do_allgather) then
1083      call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
1084 &                         mpi_enreg%comm_img,ierr)
1085    else
1086      call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
1087 &                      master_img,mpi_enreg%comm_img,ierr)
1088    end if
1089    ABI_DEALLOCATE(rbuffer)
1090    ABI_DEALLOCATE(rsize_img_all)
1091 
1092 !  Transfer buffers into gathered array_img_all (master proc only)
1093    if (use_array_all) then
1094      ABI_ALLOCATE(iimg,(mpi_enreg%nproc_img))
1095      iimg=0
1096      do jj=1,nimagetot
1097 !      The following line supposes that images are sorted by increasing index
1098        iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1
1099        ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize
1100        array_img_all(1:size1,1:size2,jj)=reshape(rbuffer_all(ibufr+1:ibufr+rsize),(/size1,size2/))
1101      end do
1102      ABI_DEALLOCATE(iimg)
1103    end if
1104 
1105 !  Free memory
1106    ABI_DEALLOCATE(rbufshft)
1107    ABI_DEALLOCATE(rbuffer_all)
1108 
1109  end if
1110 
1111 end subroutine gather_array_img_2D

m_results_img/gather_results_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  gather_results_img

FUNCTION

  Gather results_img datastructures using communicator over images (replicas) of the cell.
  Each contribution of single processor is gathered into a big array on master processor

INPUTS

  allgather= --optional, default=false--  if TRUE do ALL_GATHER instead of GATHER
  master= --optional, default=0-- index of master proc receiving gathered data (if allgather=false)
  mpi_enreg=information about MPI parallelization
  only_one_per_img= --optional, default=true--  if TRUE, the gather operation
                    is only done by one proc per image (master of the comm_cell)
  results_img(:)=<type(results_img_type)>=results_img datastructure array on each proc

SIDE EFFECTS

  results_img_all(:)=<type(results_img_type)>=global (gathered) results_img datastructure array

PARENTS

      prtimg

CHILDREN

      mkrdim,xred2xcart

SOURCE

521 subroutine gather_results_img(mpi_enreg,results_img,results_img_all,&
522 &                 master,allgather,only_one_per_img) ! optional arguments
523 
524 
525 !This section has been created automatically by the script Abilint (TD).
526 !Do not modify the following lines by hand.
527 #undef ABI_FUNC
528 #define ABI_FUNC 'gather_results_img'
529 !End of the abilint section
530 
531  implicit none
532 
533 !Arguments ------------------------------------
534 !scalars
535  integer,optional,intent(in) :: master
536  logical,optional,intent(in) :: allgather,only_one_per_img
537  type(MPI_type),intent(in) :: mpi_enreg
538 !arrays
539  type(results_img_type),intent(inout) :: results_img(:)
540  type(results_img_type),intent(inout) :: results_img_all(:)
541 
542 !Local variables-------------------------------
543 !scalars
544  integer :: ibufr,ierr,iproc,jj
545  integer :: master_all,master_img,master_one_img,natom,ngrvdw,nimage,nimagetot,npspalch,nsppol,ntypalch,ntypat
546  integer :: rsize,rsize_img
547  logical :: do_allgather,i_am_master,one_per_img,use_results_all
548  !character(len=500) :: msg
549 !arrays
550  integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:)
551  real(dp),allocatable :: rbuffer(:),rbuffer_all(:)
552 
553 !************************************************************************
554 
555  !@results_img_type
556 
557  one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img
558  do_allgather=.false.;if (present(allgather)) do_allgather=allgather
559  master_all=0;if (present(master)) master_all=master
560  i_am_master=(mpi_enreg%me==master_all)
561 
562  master_img=0;master_one_img=0
563  use_results_all= &
564 &  (((     do_allgather).and.(     one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. &
565 &   ((     do_allgather).and.(.not.one_per_img))                                            .or. &
566 &   ((.not.do_allgather).and.(     one_per_img).and.(mpi_enreg%me==master_all))             .or. &
567 &   ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img)))
568 
569 !Create global results_img_all datastructure
570  if (use_results_all) then
571    call init_results_img(results_img(1)%natom,results_img(1)%npspalch,results_img(1)%nsppol,&
572 &    results_img(1)%ntypalch,results_img(1)%ntypat,results_img_all)
573  end if
574 
575  if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then
576 
577 !  Simple copy in case of 1 image
578    if (use_results_all) then
579      if (size(results_img_all,1)<=1) then
580        call copy_results_img(results_img(1),results_img_all(1))
581        return
582      end if
583    endif
584 
585 !  Gather number of images treated by each proc
586    nimage=size(results_img,1)
587    ABI_ALLOCATE(nimage_all,(mpi_enreg%nproc_img))
588    call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr)
589    nimagetot=sum(nimage_all)
590    if (use_results_all) then
591      if (size(results_img_all,1)/=nimagetot) then
592        MSG_BUG('Wrong results_img_all size !')
593      end if
594    end if
595 
596 !  Copy natom from distributed results_img to gathered one
597    natom =results_img(1)%natom
598    npspalch =results_img(1)%npspalch
599    nsppol =results_img(1)%nsppol
600    ntypalch =results_img(1)%ntypalch
601    ntypat   =results_img(1)%ntypat
602    ngrvdw=results_img(1)%results_gs%ngrvdw
603    if (use_results_all) then
604      do jj=1,nimagetot
605        results_img_all(jj)%natom =natom
606        results_img_all(jj)%npspalch =npspalch
607        results_img_all(jj)%nsppol =nsppol
608        results_img_all(jj)%ntypalch =ntypalch
609        results_img_all(jj)%ntypat   =ntypat
610        results_img_all(jj)%results_gs%ngrvdw=ngrvdw
611      enddo
612    end if
613 
614 !  Compute number of data
615    rsize=38+n_energies+24*natom+3*nsppol+ntypat+npspalch*ntypalch+3*ngrvdw
616    rsize_img=nimage*rsize
617    ABI_ALLOCATE(rsize_img_all,(mpi_enreg%nproc_img))
618    rsize_img_all(:)=rsize*nimage_all(:)
619    ABI_DEALLOCATE(nimage_all)
620 
621 !  Compute shifts in buffer arrays for each proc
622    ABI_ALLOCATE(rbufshft,(mpi_enreg%nproc_img))
623    rbufshft(1)=0
624    do jj=2,mpi_enreg%nproc_img
625      rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1)
626    end do
627 
628 !  Load buffers
629    ABI_ALLOCATE(rbuffer,(rsize_img))
630    ibufr=0
631    do jj=1,nimage
632      rbuffer(ibufr+1)  =results_img(jj)%results_gs%deltae
633      rbuffer(ibufr+2)  =results_img(jj)%results_gs%diffor
634      rbuffer(ibufr+3)  =results_img(jj)%results_gs%entropy
635      rbuffer(ibufr+4)  =results_img(jj)%results_gs%etotal
636      rbuffer(ibufr+5)  =results_img(jj)%results_gs%fermie
637      rbuffer(ibufr+6)  =results_img(jj)%results_gs%residm
638      rbuffer(ibufr+7)  =results_img(jj)%results_gs%res2
639      rbuffer(ibufr+8)  =results_img(jj)%results_gs%vxcavg
640      rbuffer(ibufr+9:ibufr+11) =results_img(jj)%results_gs%pel(1:3)
641      rbuffer(ibufr+12:ibufr+17)=results_img(jj)%results_gs%strten(1:6)
642      rbuffer(ibufr+18:ibufr+20)=results_img(jj)%acell(1:3)
643      rbuffer(ibufr+21:ibufr+23)=results_img(jj)%rprim(1:3,1)
644      rbuffer(ibufr+24:ibufr+26)=results_img(jj)%rprim(1:3,2)
645      rbuffer(ibufr+27:ibufr+29)=results_img(jj)%rprim(1:3,3)
646      rbuffer(ibufr+30:ibufr+32)=results_img(jj)%vel_cell(1:3,1)
647      rbuffer(ibufr+33:ibufr+35)=results_img(jj)%vel_cell(1:3,2)
648      rbuffer(ibufr+36:ibufr+38)=results_img(jj)%vel_cell(1:3,3)
649      ibufr=ibufr+38
650      call energies_to_array(results_img(jj)%results_gs%energies,&
651 &                           rbuffer(ibufr+1:ibufr+n_energies),1)
652      ibufr=ibufr+n_energies
653      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%fcart(1:3,1:natom),(/3*natom/))
654      ibufr=ibufr+3*natom
655      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%fred(1:3,1:natom),(/3*natom/))
656      ibufr=ibufr+3*natom
657      rbuffer(ibufr+1:ibufr+3*nsppol)=reshape(results_img(jj)%results_gs%gaps(1:3,1:nsppol),(/3*nsppol/))
658      ibufr=ibufr+3*nsppol
659      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%gresid(1:3,1:natom),(/3*natom/))
660      ibufr=ibufr+3*natom
661      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grewtn(1:3,1:natom),(/3*natom/))
662      ibufr=ibufr+3*natom
663      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grxc(1:3,1:natom),(/3*natom/))
664      ibufr=ibufr+3*natom
665      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%synlgr(1:3,1:natom),(/3*natom/))
666      ibufr=ibufr+3*natom
667      rbuffer(ibufr+1:ibufr+ntypat)=results_img(jj)%amu(1:ntypat)
668      ibufr=ibufr+ntypat
669      rbuffer(ibufr+1:ibufr+npspalch*ntypalch)=reshape(results_img(jj)%mixalch(1:npspalch,1:ntypalch),(/npspalch*ntypalch/))
670      ibufr=ibufr+npspalch*ntypalch
671      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%xred(1:3,1:natom),(/3*natom/))
672      ibufr=ibufr+3*natom
673      rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%vel(1:3,1:natom),(/3*natom/))
674      ibufr=ibufr+3*natom
675      if (ngrvdw>0) then
676        rbuffer(ibufr+1:ibufr+3*ngrvdw)= &
677 &        reshape(results_img(jj)%results_gs%grvdw(1:3,1:ngrvdw),(/3*ngrvdw/))
678        ibufr=ibufr+3*ngrvdw
679      end if
680    end do
681    if (ibufr/=rsize_img) then
682      MSG_BUG('wrong buffer size !')
683    end if
684 
685 !  Gather all data
686    if (use_results_all)  then
687      ABI_ALLOCATE(rbuffer_all,(rsize*nimagetot))
688    end if
689    if (.not.use_results_all)  then
690      ABI_ALLOCATE(rbuffer_all,(0))
691    end if
692    if (do_allgather) then
693      call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
694 &                         mpi_enreg%comm_img,ierr)
695    else
696      call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,&
697 &                          master_img,mpi_enreg%comm_img,ierr)
698    end if
699    ABI_DEALLOCATE(rbuffer)
700    ABI_DEALLOCATE(rsize_img_all)
701 
702 !  Transfer buffers into gathered results_img_all (master proc only)
703    if (use_results_all) then
704      ABI_ALLOCATE(iimg,(mpi_enreg%nproc_img))
705      iimg=0
706      do jj=1,nimagetot
707 !      The following line supposes that images are sorted by increasing index
708        iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1
709        ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize
710        results_img_all(jj)%results_gs%deltae     =rbuffer_all(ibufr+1)
711        results_img_all(jj)%results_gs%diffor     =rbuffer_all(ibufr+2)
712        results_img_all(jj)%results_gs%entropy    =rbuffer_all(ibufr+3)
713        results_img_all(jj)%results_gs%etotal     =rbuffer_all(ibufr+4)
714        results_img_all(jj)%results_gs%fermie     =rbuffer_all(ibufr+5)
715        results_img_all(jj)%results_gs%residm     =rbuffer_all(ibufr+6)
716        results_img_all(jj)%results_gs%res2       =rbuffer_all(ibufr+7)
717        results_img_all(jj)%results_gs%vxcavg     =rbuffer_all(ibufr+8)
718        results_img_all(jj)%results_gs%pel(1:3)   =rbuffer_all(ibufr+9:ibufr+11)
719        results_img_all(jj)%results_gs%strten(1:6)=rbuffer_all(ibufr+12:ibufr+17)
720        results_img_all(jj)%acell(1:3)   =rbuffer_all(ibufr+18:ibufr+20)
721        results_img_all(jj)%rprim(1:3,1)=rbuffer_all(ibufr+21:ibufr+23)
722        results_img_all(jj)%rprim(1:3,2)=rbuffer_all(ibufr+24:ibufr+26)
723        results_img_all(jj)%rprim(1:3,3)=rbuffer_all(ibufr+27:ibufr+29)
724        results_img_all(jj)%vel_cell(1:3,1)=rbuffer_all(ibufr+30:ibufr+32)
725        results_img_all(jj)%vel_cell(1:3,2)=rbuffer_all(ibufr+33:ibufr+35)
726        results_img_all(jj)%vel_cell(1:3,3)=rbuffer_all(ibufr+36:ibufr+38)
727        ibufr=ibufr+38
728        call energies_to_array(results_img_all(jj)%results_gs%energies,&
729 &                             rbuffer_all(ibufr+1:ibufr+n_energies),-1)
730        ibufr=ibufr+n_energies
731        results_img_all(jj)%amu(1:ntypat)=rbuffer_all(ibufr+1:ibufr+ntypat)
732        results_img_all(jj)%results_gs%fcart(1:3,1:natom)= &
733 &             reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
734        ibufr=ibufr+3*natom
735        results_img_all(jj)%results_gs%fred(1:3,1:natom)= &
736 &             reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
737        ibufr=ibufr+3*natom
738        results_img_all(jj)%results_gs%gaps(1:3,1:nsppol)= &
739 &             reshape(rbuffer_all(ibufr+1:ibufr+3*nsppol),(/3,nsppol/))
740        ibufr=ibufr+3*nsppol
741        results_img_all(jj)%results_gs%gresid(1:3,1:natom)= &
742   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
743        ibufr=ibufr+3*natom
744        results_img_all(jj)%results_gs%grewtn(1:3,1:natom)= &
745   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
746        ibufr=ibufr+3*natom
747        results_img_all(jj)%results_gs%grxc(1:3,1:natom)= &
748   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
749        ibufr=ibufr+3*natom
750        results_img_all(jj)%results_gs%synlgr(1:3,1:natom)= &
751   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
752        ibufr=ibufr+3*natom
753        results_img_all(jj)%amu(1:ntypat)=rbuffer_all(ibufr+1:ibufr+ntypat)
754        ibufr=ibufr+ntypat
755        results_img_all(jj)%mixalch(1:npspalch,1:ntypalch)= &
756   &            reshape(rbuffer_all(ibufr+1:ibufr+npspalch*ntypalch),(/npspalch,ntypalch/))
757        ibufr=ibufr+npspalch*ntypalch
758        results_img_all(jj)%xred(1:3,1:natom)= &
759   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
760        ibufr=ibufr+3*natom
761        results_img_all(jj)%vel(1:3,1:natom)= &
762   &            reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/))
763        ibufr=ibufr+3*natom
764        if (ngrvdw>0) then
765          results_img_all(jj)%results_gs%grvdw(1:3,1:ngrvdw)= &
766   &              reshape(rbuffer_all(ibufr+1:ibufr+3*ngrvdw),(/3,ngrvdw/))
767          ibufr=ibufr+3*ngrvdw
768        end if
769      end do
770      ABI_DEALLOCATE(iimg)
771    end if
772 
773 !  Free memory
774    ABI_DEALLOCATE(rbufshft)
775    ABI_DEALLOCATE(rbuffer_all)
776 
777  end if
778 
779 end subroutine gather_results_img

m_results_img/get_geometry_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  get_geometry_img

FUNCTION

  From a results_img datastructure, get geometry related data

INPUTS

  natom=number of atoms
  nimage=number of images (including static ones)
  results_img(nimage)=datastructure that hold data for each image
                      (positions, forces, energy, ...)

OUTPUT

  etotal(3,natom,nimage)=total energy of each image
  fcart(3,natom,nimage)=cartesian forces in each image
  rprimd(3,3,nimage)   =dimensional primitive translations for each image
  xcart(3,natom,nimage)=cartesian coordinates of atoms in each image
  xred(3,natom,nimage) =reduced coordinates of atoms in each image

SIDE EFFECTS

PARENTS

      predict_neb,predict_steepest,predict_string

CHILDREN

      mkrdim,xred2xcart

SOURCE

1306 subroutine get_geometry_img(etotal,natom,nimage,results_img,fcart,rprimd,xcart,xred)
1307 
1308 
1309 !This section has been created automatically by the script Abilint (TD).
1310 !Do not modify the following lines by hand.
1311 #undef ABI_FUNC
1312 #define ABI_FUNC 'get_geometry_img'
1313 !End of the abilint section
1314 
1315  implicit none
1316 
1317 !Arguments ------------------------------------
1318 !scalars
1319  integer,intent(in) :: natom,nimage
1320 !arrays
1321  real(dp),intent(out) :: etotal(nimage),fcart(3,natom,nimage),rprimd(3,3,nimage)
1322  real(dp),intent(out) :: xcart(3,natom,nimage),xred(3,natom,nimage)
1323  type(results_img_type),intent(in) :: results_img(nimage)
1324 !Local variables-------------------------------
1325 !scalars
1326  integer :: iimage
1327 !arrays
1328  real(dp) :: acell(3),rprim(3,3)
1329 
1330 !************************************************************************
1331 
1332  do iimage=1,nimage
1333    acell(:)  =results_img(iimage)%acell(:)
1334    rprim(:,:)=results_img(iimage)%rprim(:,:)
1335    xred (:,:,iimage)=results_img(iimage)%xred(:,:)
1336    fcart(:,:,iimage)=results_img(iimage)%results_gs%fcart(:,:)
1337    call mkrdim(acell,rprim,rprimd(:,:,iimage))
1338    call xred2xcart(natom,rprimd(:,:,iimage),xcart(:,:,iimage),xred(:,:,iimage))
1339    etotal(iimage)=results_img(iimage)%results_gs%etotal
1340  end do
1341 
1342 end subroutine get_geometry_img

m_results_img/init_results_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  init_results_img

FUNCTION

  Init all scalars and allocatables in an array of results_img datastructures

INPUTS

  natom=number of atoms in cell
  npspalch = number of pseudopotentials for alchemical mixing  for this image
  nsppol   = number of spin channels
  ntypalch = The number of alchemical pseudoatoms  for this image

OUTPUT

SIDE EFFECTS

  results_img(:)=<type(results_img_type)>=results_img datastructure array

PARENTS

      gstateimg,m_results_img

CHILDREN

      mkrdim,xred2xcart

SOURCE

167 subroutine init_results_img(natom,npspalch,nsppol,ntypalch,ntypat,results_img)
168 
169 
170 !This section has been created automatically by the script Abilint (TD).
171 !Do not modify the following lines by hand.
172 #undef ABI_FUNC
173 #define ABI_FUNC 'init_results_img'
174 !End of the abilint section
175 
176  implicit none
177 
178 !Arguments ------------------------------------
179 !scalars
180  integer,intent(in) :: natom,npspalch,nsppol,ntypalch,ntypat
181 !arrays
182  type(results_img_type),intent(inout) :: results_img(:)
183 !Local variables-------------------------------
184 !scalars
185  integer :: ii,results_img_size
186 !arrays
187 
188 !************************************************************************
189 
190  !@results_img_type
191 
192  results_img_size=size(results_img)
193 
194  if (results_img_size>0) then
195 
196    do ii=1,results_img_size
197 
198      results_img(ii)%natom  =natom
199      results_img(ii)%npspalch  =npspalch
200      results_img(ii)%nsppol    =nsppol
201      results_img(ii)%ntypalch  =ntypalch
202      results_img(ii)%ntypat    =ntypat
203 
204      ABI_DATATYPE_ALLOCATE(results_img(ii)%results_gs,)
205      call init_results_gs(natom,nsppol,results_img(ii)%results_gs)
206 
207      ABI_ALLOCATE(results_img(ii)%acell,(3))
208      results_img(ii)%acell=zero
209      ABI_ALLOCATE(results_img(ii)%amu,(ntypat))
210      results_img(ii)%amu=zero
211      ABI_ALLOCATE(results_img(ii)%mixalch,(npspalch,ntypalch))
212      results_img(ii)%mixalch=zero
213      ABI_ALLOCATE(results_img(ii)%rprim,(3,3))
214      results_img(ii)%rprim=zero
215      ABI_ALLOCATE(results_img(ii)%xred,(3,natom))
216      results_img(ii)%xred =zero
217      ABI_ALLOCATE(results_img(ii)%vel,(3,natom))
218      results_img(ii)%vel  =zero
219      ABI_ALLOCATE(results_img(ii)%vel_cell,(3,3))
220      results_img(ii)%vel_cell=zero
221 
222    end do
223  end if
224 
225 end subroutine init_results_img

m_results_img/nullify_results_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  nullify_results_img

FUNCTION

  Nullify an array of results_img datastructures

INPUTS

OUTPUT

SIDE EFFECTS

  results_img(:)=<type(results_img_type)>=results_img datastructure array

PARENTS

CHILDREN

      mkrdim,xred2xcart

SOURCE

341 subroutine nullify_results_img(results_img)
342 
343 
344 !This section has been created automatically by the script Abilint (TD).
345 !Do not modify the following lines by hand.
346 #undef ABI_FUNC
347 #define ABI_FUNC 'nullify_results_img'
348 !End of the abilint section
349 
350  implicit none
351 
352 !Arguments ------------------------------------
353 !arrays
354  type(results_img_type),intent(inout) :: results_img(:)
355 !Local variables-------------------------------
356 !scalars
357  integer :: ii,results_img_size
358 
359 !************************************************************************
360 
361  !@results_img_type
362 
363  results_img_size=size(results_img)
364  if (results_img_size>0) then
365 
366    do ii=1,results_img_size
367      results_img(ii)%natom=0
368      results_img(ii)%npspalch=0
369      results_img(ii)%nsppol=0
370      results_img(ii)%ntypalch=0
371      results_img(ii)%ntypat=0
372      nullify(results_img(ii)%results_gs)
373    end do
374 
375  end if
376 
377 end subroutine nullify_results_img

m_results_img/results_img_type [ Types ]

[ Top ] [ m_results_img ] [ Types ]

NAME

 results_img_type

FUNCTION

 This structured datatype contains the results of a GS calculation
 for a given image of the cell:
   energy, forces, stresses,positions, velocities, cell parameter, alchemical mixing parameters

SOURCE

 77  type, public :: results_img_type
 78 
 79 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 80 ! declared in another part of ABINIT, that might need to take into account your modification.
 81 
 82 ! Integer scalar
 83 
 84   integer :: natom
 85    ! The number of atoms for this image
 86 
 87   integer :: npspalch
 88    ! The number of pseudopotentials for alchemical mixing  for this image
 89 
 90   integer :: nsppol
 91    ! The number of spin channels for this image
 92 
 93   integer :: ntypat
 94    ! The number of types of atoms for this image
 95 
 96   integer :: ntypalch
 97    ! The number of alchemical pseudoatoms  for this image
 98 
 99 ! Real (real(dp)) arrays
100 
101   real(dp), allocatable :: acell(:)
102    ! acell(3)
103    ! Dimensions of the cell
104 
105   real(dp), allocatable :: amu(:)
106    ! amu(ntypat)
107    ! Atomic masses of the different types of atoms
108 
109   real(dp), allocatable :: mixalch(:,:)
110    ! mixalch(npspalch,ntypalch)
111    ! Alchemical mixing factors
112 
113   real(dp), allocatable :: rprim(:,:)
114    ! rprim(3,3)
115    ! Primitive translations of the cell
116 
117   real(dp), allocatable :: vel(:,:)
118    ! vel(3,natom)
119    ! Velocities of the atoms
120 
121   real(dp), allocatable :: vel_cell(:,:)
122    ! vel_cell(3,3)
123    ! Time derivative of the cell parameters
124 
125   real(dp), allocatable :: xred(:,:)
126    ! xred(3,natom)
127    ! Reduced coordinates of the atoms
128 
129 ! Other types of data
130   type(results_gs_type), pointer :: results_gs => null()
131    ! Energies, forces and stresses from the GS calculation
132 
133  end type results_img_type

m_results_img/scatter_array_img [ Functions ]

[ Top ] [ m_results_img ] [ Functions ]

NAME

  scatter_array_img

FUNCTION

  Scatter an real 2D-array (part of a results_img datastructure) using communicator
  over images (replicas) of the cell.
  A big array on master processor is scattered into each contribution
  of single processor is gathered

INPUTS

  master= --optional, default=0-- index of master proc sending data
  mpi_enreg=information about MPI parallelization
  only_one_per_img= --optional, default=true--  if TRUE, the scatter operation
                    is only done by one proc per image (master of the comm_cell)
  array_img_all(:,:,:)= (real) global 2D-array (has 3 dimensions; the 3rd one is nimagetot)

SIDE EFFECTS

  array_img(:,:,:)= (real) distributed 2D-array
                    (has 3 dimensions; the 3rd one is nimage)

PARENTS

      predict_pimd

CHILDREN

      mkrdim,xred2xcart

SOURCE

1145 subroutine scatter_array_img(array_img,array_img_all,mpi_enreg,&
1146 &                            master,only_one_per_img) ! optional arguments
1147 
1148 
1149 !This section has been created automatically by the script Abilint (TD).
1150 !Do not modify the following lines by hand.
1151 #undef ABI_FUNC
1152 #define ABI_FUNC 'scatter_array_img'
1153 !End of the abilint section
1154 
1155  implicit none
1156 
1157 !Arguments ------------------------------------
1158 !scalars
1159  integer,optional,intent(in) :: master
1160  logical,optional,intent(in) :: only_one_per_img
1161  type(MPI_type),intent(in) :: mpi_enreg
1162 !arrays
1163  real(dp),intent(inout) :: array_img(:,:,:)
1164  real(dp),intent(in) :: array_img_all(:,:,:)
1165 
1166 !Local variables-------------------------------
1167 !scalars
1168  integer :: ibufr,ierr,iproc,jj
1169  integer :: master_all,master_img,master_one_img,nimage,nimagetot
1170  integer :: rsize,rsize_img,size1,size2
1171  logical :: i_am_master,one_per_img,use_array_all
1172  !character(len=500) :: msg
1173 !arrays
1174  integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:)
1175  real(dp),allocatable :: rbuffer(:),rbuffer_all(:)
1176 
1177 !************************************************************************
1178 
1179  !@results_img_type
1180 
1181  one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img
1182  master_all=0;if (present(master)) master_all=master
1183  i_am_master=(mpi_enreg%me==master_all)
1184 
1185  use_array_all=i_am_master
1186  master_img=0;master_one_img=0
1187 
1188  size1=size(array_img,1);size2=size(array_img,2)
1189  if (use_array_all) then
1190    if (size(array_img_all,1)/=size1.or.size(array_img_all,2)/=size2) then
1191      MSG_BUG('Wrong array_img_all size (1)!')
1192    end if
1193  end if
1194 
1195  if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then
1196 
1197 !  Compute (by gather operation) total number of images
1198    nimage=size(array_img,3)
1199    ABI_ALLOCATE(nimage_all,(mpi_enreg%nproc_img))
1200    call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr)
1201    nimagetot=sum(nimage_all)
1202    if (use_array_all) then
1203      if (size(array_img_all,3)/=nimagetot) then
1204        MSG_BUG('Wrong array_img_all size (2)!')
1205      endif
1206    end if
1207 
1208 !  Simple copy in case of 1 image
1209    if (nimagetot<=1) then
1210      if (use_array_all) array_img(:,:,1)=array_img_all(:,:,1)
1211 
1212    else
1213 
1214 !    Compute number of data
1215      rsize=size1*size2;rsize_img=nimage*rsize
1216      ABI_ALLOCATE(rsize_img_all,(mpi_enreg%nproc_img))
1217      rsize_img_all(:)=rsize*nimage_all(:)
1218 
1219 !    Compute shifts in buffer arrays for each proc
1220      ABI_ALLOCATE(rbufshft,(mpi_enreg%nproc_img))
1221      rbufshft(1)=0
1222      do jj=2,mpi_enreg%nproc_img
1223        rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1)
1224      end do
1225 
1226 !    Load buffer
1227      if (use_array_all)  then
1228        ABI_ALLOCATE(rbuffer_all,(rsize*nimagetot))
1229      end if
1230      if (.not.use_array_all)  then
1231        ABI_ALLOCATE(rbuffer_all,(0))
1232      end if
1233      if (use_array_all) then
1234        ABI_ALLOCATE(iimg,(mpi_enreg%nproc_img))
1235        iimg=0
1236        do jj=1,nimagetot
1237 !        The following line supposes that images are sorted by increasing index
1238          iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1
1239          ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize
1240          rbuffer_all(ibufr+1:ibufr+rsize)=reshape(array_img_all(1:size1,1:size2,jj),(/rsize/))
1241        end do
1242        ABI_DEALLOCATE(iimg)
1243       end if
1244 
1245 !    Scatter all data
1246      ABI_ALLOCATE(rbuffer,(rsize_img))
1247      call xmpi_scatterv(rbuffer_all,rsize_img_all,rbufshft,rbuffer,rsize_img,&
1248 &                       master_img,mpi_enreg%comm_img,ierr)
1249      ABI_DEALLOCATE(rbuffer_all)
1250      ABI_DEALLOCATE(rbufshft)
1251      ABI_DEALLOCATE(rsize_img_all)
1252 
1253 !    Transfered distributed buffers into array_img (master proc only)
1254      ibufr=0
1255      do jj=1,nimage
1256        array_img(1:size1,1:size2,jj)=reshape(rbuffer(ibufr+1:ibufr+rsize),(/size1,size2/))
1257        ibufr=ibufr+rsize
1258      end do
1259      ABI_DEALLOCATE(rbuffer)
1260 
1261    end if ! nimagetot<=1
1262    ABI_DEALLOCATE(nimage_all)
1263 
1264 !  Now, is requested dispatch data inside each image
1265    if (.not.one_per_img) then
1266      call xmpi_bcast(array_img,0,mpi_enreg%comm_cell,ierr)
1267    end if
1268 
1269  end if
1270 
1271 end subroutine scatter_array_img