TABLE OF CONTENTS
- ABINIT/m_results_img
- m_results_img/copy_results_img
- m_results_img/destroy_results_img
- m_results_img/gather_array_img_1D
- m_results_img/gather_array_img_2D
- m_results_img/gather_results_img
- m_results_img/get_geometry_img
- m_results_img/init_results_img
- m_results_img/nullify_results_img
- m_results_img/results_img_type
- m_results_img/scatter_array_img
ABINIT/m_results_img [ 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-2022 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
SOURCE
21 #if defined HAVE_CONFIG_H 22 #include "config.h" 23 #endif 24 25 #include "abi_common.h" 26 27 MODULE m_results_img 28 29 use defs_basis 30 use m_energies 31 use m_abicore 32 use m_results_gs 33 use m_errors 34 use m_xmpi 35 36 use defs_abitypes, only : mpi_type 37 use m_geometry, only : mkrdim, xred2xcart 38 39 implicit none 40 41 private 42 43 ! public procedures. 44 public :: init_results_img 45 public :: destroy_results_img 46 public :: nullify_results_img 47 public :: copy_results_img 48 public :: gather_results_img 49 public :: gather_array_img 50 public :: scatter_array_img 51 public :: get_geometry_img 52 53 interface gather_array_img 54 module procedure gather_array_img_1D 55 module procedure gather_array_img_2D 56 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
SOURCE
354 subroutine copy_results_img(results_img_in,results_img_out) 355 356 !Arguments ------------------------------------ 357 !arrays 358 type(results_img_type),intent(in) :: results_img_in 359 type(results_img_type),intent(inout) :: results_img_out !vz_i 360 !Local variables------------------------------- 361 !scalars 362 integer :: natom_in,natom_out,npspalch_in,npspalch_out,ntypalch_in,ntypalch_out 363 integer :: nspden_in, nspden_out, nsppol_in, nsppol_out,ntypat_in,ntypat_out 364 365 !************************************************************************ 366 367 !@results_img_type 368 369 natom_in =results_img_in%natom 370 natom_out=results_img_out%natom 371 npspalch_in =results_img_in%npspalch 372 npspalch_out=results_img_out%npspalch 373 nspden_in =results_img_in%nspden 374 nspden_out=results_img_out%nspden 375 nsppol_in =results_img_in%nsppol 376 nsppol_out=results_img_out%nsppol 377 ntypalch_in =results_img_in%ntypalch 378 ntypalch_out=results_img_out%ntypalch 379 ntypat_in =results_img_in%ntypat 380 ntypat_out=results_img_out%ntypat 381 382 if (natom_in>natom_out) then 383 if (allocated(results_img_out%xred)) then 384 ABI_FREE(results_img_out%xred) 385 end if 386 if (allocated(results_img_out%vel)) then 387 ABI_FREE(results_img_out%vel) 388 end if 389 if (allocated(results_img_out%vel_cell)) then 390 ABI_FREE(results_img_out%vel_cell) 391 end if 392 393 if (allocated(results_img_in%xred)) then 394 ABI_MALLOC(results_img_out%xred,(3,natom_in)) 395 end if 396 if (allocated(results_img_in%vel)) then 397 ABI_MALLOC(results_img_out%vel,(3,natom_in)) 398 end if 399 if (allocated(results_img_in%vel_cell)) then 400 ABI_MALLOC(results_img_out%vel_cell,(3,3)) 401 end if 402 end if 403 404 if (npspalch_in>npspalch_out .or. ntypalch_in>ntypalch_out) then 405 if (allocated(results_img_out%mixalch)) then 406 ABI_FREE(results_img_out%mixalch) 407 end if 408 endif 409 410 if (ntypat_in>ntypat_out) then 411 if (allocated(results_img_in%amu)) then 412 ABI_MALLOC(results_img_out%amu,(ntypat_in)) 413 endif 414 endif 415 416 results_img_out%natom =results_img_in%natom 417 results_img_out%npspalch =results_img_in%npspalch 418 results_img_out%nspden =results_img_in%nspden 419 results_img_out%nsppol =results_img_in%nsppol 420 results_img_out%ntypalch =results_img_in%ntypalch 421 results_img_out%ntypat =results_img_in%ntypat 422 423 call copy_results_gs(results_img_in%results_gs,results_img_out%results_gs) 424 425 if (allocated(results_img_in%acell)) results_img_out%acell(:)=results_img_in%acell(:) 426 if (allocated(results_img_in%amu)) results_img_out%amu(1:ntypat_in)=results_img_in%amu(1:ntypat_in) 427 if (allocated(results_img_in%mixalch)) & 428 & results_img_out%mixalch(1:npspalch_in,1:ntypalch_in)=results_img_in%mixalch(1:npspalch_in,1:ntypalch_in) 429 if (allocated(results_img_in%rprim)) results_img_out%rprim(:,:)=results_img_in%rprim(:,:) 430 if (allocated(results_img_in%xred)) results_img_out%xred(:,1:natom_in)=results_img_in%xred(:,1:natom_in) 431 if (allocated(results_img_in%vel)) results_img_out%vel(:,1:natom_in)=results_img_in%vel(:,1:natom_in) 432 if (allocated(results_img_in%vel_cell))results_img_out%vel_cell(:,:)=results_img_in%vel_cell(:,:) 433 434 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
SOURCE
230 subroutine destroy_results_img(results_img) 231 232 !Arguments ------------------------------------ 233 !arrays 234 type(results_img_type),intent(inout) :: results_img(:) 235 !Local variables------------------------------- 236 !scalars 237 integer :: ii,results_img_size 238 239 !************************************************************************ 240 241 !@results_img_type 242 243 results_img_size=size(results_img) 244 if (results_img_size>0) then 245 246 do ii=1,results_img_size 247 248 results_img(ii)%natom=0 249 results_img(ii)%npspalch=0 250 results_img(ii)%nspden =0 251 results_img(ii)%nsppol =0 252 results_img(ii)%ntypalch=0 253 results_img(ii)%ntypat=0 254 255 if (associated(results_img(ii)%results_gs)) then 256 call destroy_results_gs(results_img(ii)%results_gs) 257 ABI_FREE(results_img(ii)%results_gs) 258 end if 259 260 if (allocated(results_img(ii)%acell)) then 261 ABI_FREE(results_img(ii)%acell) 262 end if 263 if (allocated(results_img(ii)%amu)) then 264 ABI_FREE(results_img(ii)%amu) 265 end if 266 if (allocated(results_img(ii)%mixalch)) then 267 ABI_FREE(results_img(ii)%mixalch) 268 end if 269 if (allocated(results_img(ii)%rprim)) then 270 ABI_FREE(results_img(ii)%rprim) 271 end if 272 if (allocated(results_img(ii)%xred)) then 273 ABI_FREE(results_img(ii)%xred) 274 end if 275 if (allocated(results_img(ii)%vel)) then 276 ABI_FREE(results_img(ii)%vel) 277 end if 278 if (allocated(results_img(ii)%vel_cell)) then 279 ABI_FREE(results_img(ii)%vel_cell) 280 end if 281 end do 282 283 end if 284 285 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)
SOURCE
755 subroutine gather_array_img_1D(array_img,array_img_all,mpi_enreg,& 756 & master,allgather,only_one_per_img) ! optional arguments 757 758 !Arguments ------------------------------------ 759 !scalars 760 integer,optional,intent(in) :: master 761 logical,optional,intent(in) :: allgather,only_one_per_img 762 type(MPI_type),intent(in) :: mpi_enreg 763 !arrays 764 real(dp),intent(in) :: array_img(:,:) 765 real(dp),intent(inout) :: array_img_all(:,:) 766 767 !Local variables------------------------------- 768 !scalars 769 integer :: ibufr,ierr,iproc,jj 770 integer :: master_all,master_img,master_one_img,nimage,nimagetot 771 integer :: rsize,rsize_img,size1 772 logical :: do_allgather,i_am_master,one_per_img,use_array_all 773 !character(len=500) :: msg 774 !arrays 775 integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:) 776 real(dp),allocatable :: rbuffer(:),rbuffer_all(:) 777 778 !************************************************************************ 779 780 !@results_img_type 781 782 one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img 783 do_allgather=.false.;if (present(allgather)) do_allgather=allgather 784 master_all=0;if (present(master)) master_all=master 785 i_am_master=(mpi_enreg%me==master_all) 786 787 master_img=0;master_one_img=0 788 use_array_all= & 789 & ((( do_allgather).and.( one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. & 790 & (( do_allgather).and.(.not.one_per_img)) .or. & 791 & ((.not.do_allgather).and.( one_per_img).and.(mpi_enreg%me==master_all)) .or. & 792 & ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img))) 793 794 size1=size(array_img,1) 795 if (use_array_all) then 796 if (size(array_img_all,1)/=size1) then 797 ABI_BUG('Wrong array_img_all size (1)') 798 end if 799 end if 800 801 if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then 802 803 ! Simple copy in case of 1 image 804 if (use_array_all) then 805 if (size(array_img_all,2)<=1) then 806 array_img_all(:,1)=array_img(:,1) 807 return 808 end if 809 endif 810 811 ! Gather number of images treated by each proc 812 nimage=size(array_img,2) 813 ABI_MALLOC(nimage_all,(mpi_enreg%nproc_img)) 814 call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr) 815 nimagetot=sum(nimage_all) 816 if (use_array_all) then 817 if (size(array_img_all,2)/=nimagetot) then 818 ABI_BUG('Wrong array_img_all size (2)!') 819 endif 820 end if 821 822 ! Compute number of data 823 rsize=size1;rsize_img=nimage*rsize 824 ABI_MALLOC(rsize_img_all,(mpi_enreg%nproc_img)) 825 rsize_img_all(:)=rsize*nimage_all(:) 826 ABI_FREE(nimage_all) 827 828 ! Compute shifts in buffer arrays for each proc 829 ABI_MALLOC(rbufshft,(mpi_enreg%nproc_img)) 830 rbufshft(1)=0 831 do jj=2,mpi_enreg%nproc_img 832 rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1) 833 end do 834 835 ! Load buffers 836 ABI_MALLOC(rbuffer,(rsize_img)) 837 ibufr=0 838 do jj=1,nimage 839 rbuffer(ibufr+1:ibufr+rsize)=reshape(array_img(1:size1,jj),(/rsize/)) 840 ibufr=ibufr+rsize 841 end do 842 843 ! Gather all data 844 if (use_array_all) then 845 ABI_MALLOC(rbuffer_all,(rsize*nimagetot)) 846 end if 847 if (.not.use_array_all) then 848 ABI_MALLOC(rbuffer_all,(0)) 849 end if 850 if (do_allgather) then 851 call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 852 & mpi_enreg%comm_img,ierr) 853 else 854 call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 855 & master_img,mpi_enreg%comm_img,ierr) 856 end if 857 ABI_FREE(rbuffer) 858 ABI_FREE(rsize_img_all) 859 860 ! Transfer buffers into gathered array_img_all (master proc only) 861 if (use_array_all) then 862 ABI_MALLOC(iimg,(mpi_enreg%nproc_img)) 863 iimg=0 864 do jj=1,nimagetot 865 ! The following line supposes that images are sorted by increasing index 866 iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1 867 ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize 868 array_img_all(1:size1,jj)=reshape(rbuffer_all(ibufr+1:ibufr+rsize),(/size1/)) 869 end do 870 ABI_FREE(iimg) 871 end if 872 873 ! Free memory 874 ABI_FREE(rbufshft) 875 ABI_FREE(rbuffer_all) 876 877 end if 878 879 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)
SOURCE
907 subroutine gather_array_img_2D(array_img,array_img_all,mpi_enreg,& 908 & master,allgather,only_one_per_img) ! optional arguments 909 910 !Arguments ------------------------------------ 911 !scalars 912 integer,optional,intent(in) :: master 913 logical,optional,intent(in) :: allgather,only_one_per_img 914 type(MPI_type),intent(in) :: mpi_enreg 915 !arrays 916 real(dp),intent(in) :: array_img(:,:,:) 917 real(dp),intent(inout) :: array_img_all(:,:,:) 918 919 !Local variables------------------------------- 920 !scalars 921 integer :: ibufr,ierr,iproc,jj 922 integer :: master_all,master_img,master_one_img,nimage,nimagetot 923 integer :: rsize,rsize_img,size1,size2 924 logical :: do_allgather,i_am_master,one_per_img,use_array_all 925 !character(len=500) :: msg 926 !arrays 927 integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:) 928 real(dp),allocatable :: rbuffer(:),rbuffer_all(:) 929 930 !************************************************************************ 931 932 !@results_img_type 933 934 one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img 935 do_allgather=.false.;if (present(allgather)) do_allgather=allgather 936 master_all=0;if (present(master)) master_all=master 937 i_am_master=(mpi_enreg%me==master_all) 938 939 master_img=0;master_one_img=0 940 use_array_all= & 941 & ((( do_allgather).and.( one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. & 942 & (( do_allgather).and.(.not.one_per_img)) .or. & 943 & ((.not.do_allgather).and.( one_per_img).and.(mpi_enreg%me==master_all)) .or. & 944 & ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img))) 945 946 size1=size(array_img,1);size2=size(array_img,2) 947 if (use_array_all) then 948 if (size(array_img_all,1)/=size1.or.size(array_img_all,2)/=size2) then 949 ABI_BUG('Wrong array_img_all size (1)!') 950 end if 951 end if 952 953 if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then 954 955 ! Simple copy in case of 1 image 956 if (use_array_all) then 957 if (size(array_img_all,3)<=1) then 958 array_img_all(:,:,1)=array_img(:,:,1) 959 return 960 end if 961 endif 962 963 ! Gather number of images treated by each proc 964 nimage=size(array_img,3) 965 ABI_MALLOC(nimage_all,(mpi_enreg%nproc_img)) 966 call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr) 967 nimagetot=sum(nimage_all) 968 if (use_array_all) then 969 if (size(array_img_all,3)/=nimagetot) then 970 ABI_BUG('Wrong array_img_all size (2)!') 971 endif 972 end if 973 974 ! Compute number of data 975 rsize=size1*size2;rsize_img=nimage*rsize 976 ABI_MALLOC(rsize_img_all,(mpi_enreg%nproc_img)) 977 rsize_img_all(:)=rsize*nimage_all(:) 978 ABI_FREE(nimage_all) 979 980 ! Compute shifts in buffer arrays for each proc 981 ABI_MALLOC(rbufshft,(mpi_enreg%nproc_img)) 982 rbufshft(1)=0 983 do jj=2,mpi_enreg%nproc_img 984 rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1) 985 end do 986 987 ! Load buffers 988 ABI_MALLOC(rbuffer,(rsize_img)) 989 ibufr=0 990 do jj=1,nimage 991 rbuffer(ibufr+1:ibufr+rsize)=reshape(array_img(1:size1,1:size2,jj),(/rsize/)) 992 ibufr=ibufr+rsize 993 end do 994 995 ! Gather all data 996 if (use_array_all) then 997 ABI_MALLOC(rbuffer_all,(rsize*nimagetot)) 998 end if 999 if (.not.use_array_all) then 1000 ABI_MALLOC(rbuffer_all,(0)) 1001 end if 1002 if (do_allgather) then 1003 call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 1004 & mpi_enreg%comm_img,ierr) 1005 else 1006 call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 1007 & master_img,mpi_enreg%comm_img,ierr) 1008 end if 1009 ABI_FREE(rbuffer) 1010 ABI_FREE(rsize_img_all) 1011 1012 ! Transfer buffers into gathered array_img_all (master proc only) 1013 if (use_array_all) then 1014 ABI_MALLOC(iimg,(mpi_enreg%nproc_img)) 1015 iimg=0 1016 do jj=1,nimagetot 1017 ! The following line supposes that images are sorted by increasing index 1018 iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1 1019 ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize 1020 array_img_all(1:size1,1:size2,jj)=reshape(rbuffer_all(ibufr+1:ibufr+rsize),(/size1,size2/)) 1021 end do 1022 ABI_FREE(iimg) 1023 end if 1024 1025 ! Free memory 1026 ABI_FREE(rbufshft) 1027 ABI_FREE(rbuffer_all) 1028 1029 end if 1030 1031 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
SOURCE
460 subroutine gather_results_img(mpi_enreg,results_img,results_img_all,& 461 & master,allgather,only_one_per_img) ! optional arguments 462 463 !Arguments ------------------------------------ 464 !scalars 465 integer,optional,intent(in) :: master 466 logical,optional,intent(in) :: allgather,only_one_per_img 467 type(MPI_type),intent(in) :: mpi_enreg 468 !arrays 469 type(results_img_type),intent(inout) :: results_img(:) 470 type(results_img_type),intent(inout) :: results_img_all(:) 471 472 !Local variables------------------------------- 473 !scalars 474 integer :: ibufr,ierr,iproc,jj 475 integer :: master_all,master_img,master_one_img 476 integer :: natom,ngrvdw,nimage,nimagetot,npspalch,nspden,nsppol,ntypalch,ntypat 477 integer :: rsize,rsize_img 478 logical :: do_allgather,i_am_master,one_per_img,use_results_all 479 !character(len=500) :: msg 480 !arrays 481 integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:) 482 real(dp),allocatable :: rbuffer(:),rbuffer_all(:) 483 484 !************************************************************************ 485 486 !@results_img_type 487 488 one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img 489 do_allgather=.false.;if (present(allgather)) do_allgather=allgather 490 master_all=0;if (present(master)) master_all=master 491 i_am_master=(mpi_enreg%me==master_all) 492 493 master_img=0;master_one_img=0 494 use_results_all= & 495 & ((( do_allgather).and.( one_per_img).and.(mpi_enreg%me_cell==master_one_img)) .or. & 496 & (( do_allgather).and.(.not.one_per_img)) .or. & 497 & ((.not.do_allgather).and.( one_per_img).and.(mpi_enreg%me==master_all)) .or. & 498 & ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_enreg%me_img==master_img))) 499 500 !Create global results_img_all datastructure 501 if (use_results_all) then 502 call init_results_img(results_img(1)%natom,results_img(1)%npspalch,results_img(1)%nspden,results_img(1)%nsppol,& 503 & results_img(1)%ntypalch,results_img(1)%ntypat,results_img_all) 504 end if 505 506 if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then 507 508 ! Simple copy in case of 1 image 509 if (use_results_all) then 510 if (size(results_img_all,1)<=1) then 511 call copy_results_img(results_img(1),results_img_all(1)) 512 return 513 end if 514 endif 515 516 ! Gather number of images treated by each proc 517 nimage=size(results_img,1) 518 ABI_MALLOC(nimage_all,(mpi_enreg%nproc_img)) 519 call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr) 520 nimagetot=sum(nimage_all) 521 if (use_results_all) then 522 if (size(results_img_all,1)/=nimagetot) then 523 ABI_BUG('Wrong results_img_all size !') 524 end if 525 end if 526 527 ! Copy natom from distributed results_img to gathered one 528 natom =results_img(1)%natom 529 npspalch =results_img(1)%npspalch 530 nspden =results_img(1)%nspden 531 nsppol =results_img(1)%nsppol 532 ntypalch =results_img(1)%ntypalch 533 ntypat =results_img(1)%ntypat 534 ngrvdw=results_img(1)%results_gs%ngrvdw 535 if (use_results_all) then 536 do jj=1,nimagetot 537 results_img_all(jj)%natom =natom 538 results_img_all(jj)%npspalch =npspalch 539 results_img_all(jj)%nspden =nspden 540 results_img_all(jj)%nsppol =nsppol 541 results_img_all(jj)%ntypalch =ntypalch 542 results_img_all(jj)%ntypat =ntypat 543 results_img_all(jj)%results_gs%ngrvdw=ngrvdw 544 enddo 545 end if 546 547 ! Compute number of data 548 rsize=38+n_energies+(30+nspden)*natom+3*nsppol+ntypat+npspalch*ntypalch+3*ngrvdw 549 rsize_img=nimage*rsize 550 ABI_MALLOC(rsize_img_all,(mpi_enreg%nproc_img)) 551 rsize_img_all(:)=rsize*nimage_all(:) 552 ABI_FREE(nimage_all) 553 554 ! Compute shifts in buffer arrays for each proc 555 ABI_MALLOC(rbufshft,(mpi_enreg%nproc_img)) 556 rbufshft(1)=0 557 do jj=2,mpi_enreg%nproc_img 558 rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1) 559 end do 560 561 ! Load buffers 562 ABI_MALLOC(rbuffer,(rsize_img)) 563 ibufr=0 564 do jj=1,nimage 565 rbuffer(ibufr+1) =results_img(jj)%results_gs%deltae 566 rbuffer(ibufr+2) =results_img(jj)%results_gs%diffor 567 rbuffer(ibufr+3) =results_img(jj)%results_gs%entropy 568 rbuffer(ibufr+4) =results_img(jj)%results_gs%etotal 569 rbuffer(ibufr+5) =results_img(jj)%results_gs%fermie 570 rbuffer(ibufr+6) =results_img(jj)%results_gs%residm 571 rbuffer(ibufr+7) =results_img(jj)%results_gs%res2 572 rbuffer(ibufr+8) =results_img(jj)%results_gs%vxcavg 573 rbuffer(ibufr+9:ibufr+11) =results_img(jj)%results_gs%pel(1:3) 574 rbuffer(ibufr+12:ibufr+17)=results_img(jj)%results_gs%strten(1:6) 575 rbuffer(ibufr+18:ibufr+20)=results_img(jj)%acell(1:3) 576 rbuffer(ibufr+21:ibufr+23)=results_img(jj)%rprim(1:3,1) 577 rbuffer(ibufr+24:ibufr+26)=results_img(jj)%rprim(1:3,2) 578 rbuffer(ibufr+27:ibufr+29)=results_img(jj)%rprim(1:3,3) 579 rbuffer(ibufr+30:ibufr+32)=results_img(jj)%vel_cell(1:3,1) 580 rbuffer(ibufr+33:ibufr+35)=results_img(jj)%vel_cell(1:3,2) 581 rbuffer(ibufr+36:ibufr+38)=results_img(jj)%vel_cell(1:3,3) 582 ibufr=ibufr+38 583 call energies_to_array(results_img(jj)%results_gs%energies,& 584 & rbuffer(ibufr+1:ibufr+n_energies),1) 585 ibufr=ibufr+n_energies 586 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%fcart(1:3,1:natom),(/3*natom/)) 587 ibufr=ibufr+3*natom 588 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%gred(1:3,1:natom),(/3*natom/)) 589 ibufr=ibufr+3*natom 590 rbuffer(ibufr+1:ibufr+3*nsppol)=reshape(results_img(jj)%results_gs%gaps(1:3,1:nsppol),(/3*nsppol/)) 591 ibufr=ibufr+3*nsppol 592 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grchempottn(1:3,1:natom),(/3*natom/)) 593 ibufr=ibufr+3*natom 594 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grcondft(1:3,1:natom),(/3*natom/)) 595 ibufr=ibufr+3*natom 596 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%gresid(1:3,1:natom),(/3*natom/)) 597 ibufr=ibufr+3*natom 598 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grewtn(1:3,1:natom),(/3*natom/)) 599 ibufr=ibufr+3*natom 600 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%grxc(1:3,1:natom),(/3*natom/)) 601 ibufr=ibufr+3*natom 602 rbuffer(ibufr+1:ibufr+nspden*natom)=reshape(results_img(jj)%results_gs%intgres(1:nspden,1:natom),(/nspden*natom/)) 603 ibufr=ibufr+nspden*natom 604 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%results_gs%synlgr(1:3,1:natom),(/3*natom/)) 605 ibufr=ibufr+3*natom 606 rbuffer(ibufr+1:ibufr+ntypat)=results_img(jj)%amu(1:ntypat) 607 ibufr=ibufr+ntypat 608 rbuffer(ibufr+1:ibufr+npspalch*ntypalch)=reshape(results_img(jj)%mixalch(1:npspalch,1:ntypalch),(/npspalch*ntypalch/)) 609 ibufr=ibufr+npspalch*ntypalch 610 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%xred(1:3,1:natom),(/3*natom/)) 611 ibufr=ibufr+3*natom 612 rbuffer(ibufr+1:ibufr+3*natom)=reshape(results_img(jj)%vel(1:3,1:natom),(/3*natom/)) 613 ibufr=ibufr+3*natom 614 if (ngrvdw>0) then 615 rbuffer(ibufr+1:ibufr+3*ngrvdw)= & 616 & reshape(results_img(jj)%results_gs%grvdw(1:3,1:ngrvdw),(/3*ngrvdw/)) 617 ibufr=ibufr+3*ngrvdw 618 end if 619 end do 620 if (ibufr/=rsize_img) then 621 ABI_BUG('wrong buffer size !') 622 end if 623 624 ! Gather all data 625 if (use_results_all) then 626 ABI_MALLOC(rbuffer_all,(rsize*nimagetot)) 627 end if 628 if (.not.use_results_all) then 629 ABI_MALLOC(rbuffer_all,(0)) 630 end if 631 if (do_allgather) then 632 call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 633 & mpi_enreg%comm_img,ierr) 634 else 635 call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 636 & master_img,mpi_enreg%comm_img,ierr) 637 end if 638 ABI_FREE(rbuffer) 639 ABI_FREE(rsize_img_all) 640 641 ! Transfer buffers into gathered results_img_all (master proc only) 642 if (use_results_all) then 643 ABI_MALLOC(iimg,(mpi_enreg%nproc_img)) 644 iimg=0 645 do jj=1,nimagetot 646 ! The following line supposes that images are sorted by increasing index 647 iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1 648 ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize 649 results_img_all(jj)%results_gs%deltae =rbuffer_all(ibufr+1) 650 results_img_all(jj)%results_gs%diffor =rbuffer_all(ibufr+2) 651 results_img_all(jj)%results_gs%entropy =rbuffer_all(ibufr+3) 652 results_img_all(jj)%results_gs%etotal =rbuffer_all(ibufr+4) 653 results_img_all(jj)%results_gs%fermie =rbuffer_all(ibufr+5) 654 results_img_all(jj)%results_gs%residm =rbuffer_all(ibufr+6) 655 results_img_all(jj)%results_gs%res2 =rbuffer_all(ibufr+7) 656 results_img_all(jj)%results_gs%vxcavg =rbuffer_all(ibufr+8) 657 results_img_all(jj)%results_gs%pel(1:3) =rbuffer_all(ibufr+9:ibufr+11) 658 results_img_all(jj)%results_gs%strten(1:6)=rbuffer_all(ibufr+12:ibufr+17) 659 results_img_all(jj)%acell(1:3) =rbuffer_all(ibufr+18:ibufr+20) 660 results_img_all(jj)%rprim(1:3,1)=rbuffer_all(ibufr+21:ibufr+23) 661 results_img_all(jj)%rprim(1:3,2)=rbuffer_all(ibufr+24:ibufr+26) 662 results_img_all(jj)%rprim(1:3,3)=rbuffer_all(ibufr+27:ibufr+29) 663 results_img_all(jj)%vel_cell(1:3,1)=rbuffer_all(ibufr+30:ibufr+32) 664 results_img_all(jj)%vel_cell(1:3,2)=rbuffer_all(ibufr+33:ibufr+35) 665 results_img_all(jj)%vel_cell(1:3,3)=rbuffer_all(ibufr+36:ibufr+38) 666 ibufr=ibufr+38 667 call energies_to_array(results_img_all(jj)%results_gs%energies,& 668 & rbuffer_all(ibufr+1:ibufr+n_energies),-1) 669 ibufr=ibufr+n_energies 670 results_img_all(jj)%amu(1:ntypat)=rbuffer_all(ibufr+1:ibufr+ntypat) 671 results_img_all(jj)%results_gs%fcart(1:3,1:natom)= & 672 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 673 ibufr=ibufr+3*natom 674 results_img_all(jj)%results_gs%gred(1:3,1:natom)= & 675 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 676 ibufr=ibufr+3*natom 677 results_img_all(jj)%results_gs%gaps(1:3,1:nsppol)= & 678 & reshape(rbuffer_all(ibufr+1:ibufr+3*nsppol),(/3,nsppol/)) 679 ibufr=ibufr+3*nsppol 680 results_img_all(jj)%results_gs%grchempottn(1:3,1:natom)= & 681 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 682 ibufr=ibufr+3*natom 683 results_img_all(jj)%results_gs%grcondft(1:3,1:natom)= & 684 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 685 ibufr=ibufr+3*natom 686 results_img_all(jj)%results_gs%gresid(1:3,1:natom)= & 687 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 688 ibufr=ibufr+3*natom 689 results_img_all(jj)%results_gs%grewtn(1:3,1:natom)= & 690 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 691 ibufr=ibufr+3*natom 692 results_img_all(jj)%results_gs%grxc(1:3,1:natom)= & 693 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 694 ibufr=ibufr+3*natom 695 results_img_all(jj)%results_gs%intgres(1:nspden,1:natom)= & 696 & reshape(rbuffer_all(ibufr+1:ibufr+nspden*natom),(/nspden,natom/)) 697 ibufr=ibufr+nspden*natom 698 results_img_all(jj)%results_gs%synlgr(1:3,1:natom)= & 699 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 700 ibufr=ibufr+3*natom 701 results_img_all(jj)%amu(1:ntypat)=rbuffer_all(ibufr+1:ibufr+ntypat) 702 ibufr=ibufr+ntypat 703 results_img_all(jj)%mixalch(1:npspalch,1:ntypalch)= & 704 & reshape(rbuffer_all(ibufr+1:ibufr+npspalch*ntypalch),(/npspalch,ntypalch/)) 705 ibufr=ibufr+npspalch*ntypalch 706 results_img_all(jj)%xred(1:3,1:natom)= & 707 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 708 ibufr=ibufr+3*natom 709 results_img_all(jj)%vel(1:3,1:natom)= & 710 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom),(/3,natom/)) 711 ibufr=ibufr+3*natom 712 if (ngrvdw>0) then 713 results_img_all(jj)%results_gs%grvdw(1:3,1:ngrvdw)= & 714 & reshape(rbuffer_all(ibufr+1:ibufr+3*ngrvdw),(/3,ngrvdw/)) 715 ibufr=ibufr+3*ngrvdw 716 end if 717 end do 718 ABI_FREE(iimg) 719 end if 720 721 ! Free memory 722 ABI_FREE(rbufshft) 723 ABI_FREE(rbuffer_all) 724 725 end if 726 727 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
SOURCE
1205 subroutine get_geometry_img(etotal,natom,nimage,results_img,fcart,rprimd,xcart,xred) 1206 1207 !Arguments ------------------------------------ 1208 !scalars 1209 integer,intent(in) :: natom,nimage 1210 !arrays 1211 real(dp),intent(out) :: etotal(nimage),fcart(3,natom,nimage),rprimd(3,3,nimage) 1212 real(dp),intent(out) :: xcart(3,natom,nimage),xred(3,natom,nimage) 1213 type(results_img_type),intent(in) :: results_img(nimage) 1214 !Local variables------------------------------- 1215 !scalars 1216 integer :: iimage 1217 !arrays 1218 real(dp) :: acell(3),rprim(3,3) 1219 1220 !************************************************************************ 1221 1222 do iimage=1,nimage 1223 acell(:) =results_img(iimage)%acell(:) 1224 rprim(:,:)=results_img(iimage)%rprim(:,:) 1225 xred (:,:,iimage)=results_img(iimage)%xred(:,:) 1226 fcart(:,:,iimage)=results_img(iimage)%results_gs%fcart(:,:) 1227 call mkrdim(acell,rprim,rprimd(:,:,iimage)) 1228 call xred2xcart(natom,rprimd(:,:,iimage),xcart(:,:,iimage),xred(:,:,iimage)) 1229 etotal(iimage)=results_img(iimage)%results_gs%etotal 1230 end do 1231 1232 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 nspden = number of spin-density channels 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
SOURCE
159 subroutine init_results_img(natom,npspalch,nspden,nsppol,ntypalch,ntypat,results_img) 160 161 !Arguments ------------------------------------ 162 !scalars 163 integer,intent(in) :: natom,npspalch,nspden,nsppol,ntypalch,ntypat 164 !arrays 165 type(results_img_type),intent(inout) :: results_img(:) 166 !Local variables------------------------------- 167 !scalars 168 integer :: ii,results_img_size 169 !arrays 170 171 !************************************************************************ 172 173 !@results_img_type 174 175 results_img_size=size(results_img) 176 177 if (results_img_size>0) then 178 179 do ii=1,results_img_size 180 181 results_img(ii)%natom =natom 182 results_img(ii)%npspalch =npspalch 183 results_img(ii)%nspden =nspden 184 results_img(ii)%nsppol =nsppol 185 results_img(ii)%ntypalch =ntypalch 186 results_img(ii)%ntypat =ntypat 187 188 ABI_MALLOC(results_img(ii)%results_gs,) 189 call init_results_gs(natom,nspden,nsppol,results_img(ii)%results_gs) 190 191 ABI_MALLOC(results_img(ii)%acell,(3)) 192 results_img(ii)%acell=zero 193 ABI_MALLOC(results_img(ii)%amu,(ntypat)) 194 results_img(ii)%amu=zero 195 ABI_MALLOC(results_img(ii)%mixalch,(npspalch,ntypalch)) 196 results_img(ii)%mixalch=zero 197 ABI_MALLOC(results_img(ii)%rprim,(3,3)) 198 results_img(ii)%rprim=zero 199 ABI_MALLOC(results_img(ii)%xred,(3,natom)) 200 results_img(ii)%xred =zero 201 ABI_MALLOC(results_img(ii)%vel,(3,natom)) 202 results_img(ii)%vel =zero 203 ABI_MALLOC(results_img(ii)%vel_cell,(3,3)) 204 results_img(ii)%vel_cell=zero 205 206 end do 207 end if 208 209 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
SOURCE
306 subroutine nullify_results_img(results_img) 307 308 !Arguments ------------------------------------ 309 !arrays 310 type(results_img_type),intent(inout) :: results_img(:) 311 !Local variables------------------------------- 312 !scalars 313 integer :: ii,results_img_size 314 315 !************************************************************************ 316 317 !@results_img_type 318 319 results_img_size=size(results_img) 320 if (results_img_size>0) then 321 322 do ii=1,results_img_size 323 results_img(ii)%natom=0 324 results_img(ii)%npspalch=0 325 results_img(ii)%nspden=0 326 results_img(ii)%nsppol=0 327 results_img(ii)%ntypalch=0 328 results_img(ii)%ntypat=0 329 nullify(results_img(ii)%results_gs) 330 end do 331 332 end if 333 334 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
71 type, public :: results_img_type 72 73 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 74 ! declared in another part of ABINIT, that might need to take into account your modification. 75 76 ! Integer scalar 77 78 integer :: natom 79 ! The number of atoms for this image 80 81 integer :: npspalch 82 ! The number of pseudopotentials for alchemical mixing for this image 83 84 integer :: nspden 85 ! The number of spin-density channels for this image 86 87 integer :: nsppol 88 ! The number of spin channels for this image 89 90 integer :: ntypat 91 ! The number of types of atoms for this image 92 93 integer :: ntypalch 94 ! The number of alchemical pseudoatoms for this image 95 96 ! Real (real(dp)) arrays 97 98 real(dp), allocatable :: acell(:) 99 ! acell(3) 100 ! Dimensions of the cell 101 102 real(dp), allocatable :: amu(:) 103 ! amu(ntypat) 104 ! Atomic masses of the different types of atoms 105 106 real(dp), allocatable :: mixalch(:,:) 107 ! mixalch(npspalch,ntypalch) 108 ! Alchemical mixing factors 109 110 real(dp), allocatable :: rprim(:,:) 111 ! rprim(3,3) 112 ! Primitive translations of the cell 113 114 real(dp), allocatable :: vel(:,:) 115 ! vel(3,natom) 116 ! Velocities of the atoms 117 118 real(dp), allocatable :: vel_cell(:,:) 119 ! vel_cell(3,3) 120 ! Time derivative of the cell parameters 121 122 real(dp), allocatable :: xred(:,:) 123 ! xred(3,natom) 124 ! Reduced coordinates of the atoms 125 126 ! Other types of data 127 type(results_gs_type), pointer :: results_gs => null() 128 ! Energies, forces and stresses from the GS calculation 129 130 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)
SOURCE
1059 subroutine scatter_array_img(array_img,array_img_all,mpi_enreg,& 1060 & master,only_one_per_img) ! optional arguments 1061 1062 !Arguments ------------------------------------ 1063 !scalars 1064 integer,optional,intent(in) :: master 1065 logical,optional,intent(in) :: only_one_per_img 1066 type(MPI_type),intent(in) :: mpi_enreg 1067 !arrays 1068 real(dp),intent(inout) :: array_img(:,:,:) 1069 real(dp),intent(in) :: array_img_all(:,:,:) 1070 1071 !Local variables------------------------------- 1072 !scalars 1073 integer :: ibufr,ierr,iproc,jj 1074 integer :: master_all,master_img,master_one_img,nimage,nimagetot 1075 integer :: rsize,rsize_img,size1,size2 1076 logical :: i_am_master,one_per_img,use_array_all 1077 !character(len=500) :: msg 1078 !arrays 1079 integer,allocatable :: iimg(:),nimage_all(:),rbufshft(:),rsize_img_all(:) 1080 real(dp),allocatable :: rbuffer(:),rbuffer_all(:) 1081 1082 !************************************************************************ 1083 1084 !@results_img_type 1085 1086 one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img 1087 master_all=0;if (present(master)) master_all=master 1088 i_am_master=(mpi_enreg%me==master_all) 1089 1090 use_array_all=i_am_master 1091 master_img=0;master_one_img=0 1092 1093 size1=size(array_img,1);size2=size(array_img,2) 1094 if (use_array_all) then 1095 if (size(array_img_all,1)/=size1.or.size(array_img_all,2)/=size2) then 1096 ABI_BUG('Wrong array_img_all size (1)!') 1097 end if 1098 end if 1099 1100 if ((.not.one_per_img).or.(mpi_enreg%me_cell==master_one_img)) then 1101 1102 ! Compute (by gather operation) total number of images 1103 nimage=size(array_img,3) 1104 ABI_MALLOC(nimage_all,(mpi_enreg%nproc_img)) 1105 call xmpi_allgather(nimage,nimage_all,mpi_enreg%comm_img,ierr) 1106 nimagetot=sum(nimage_all) 1107 if (use_array_all) then 1108 if (size(array_img_all,3)/=nimagetot) then 1109 ABI_BUG('Wrong array_img_all size (2)!') 1110 endif 1111 end if 1112 1113 ! Simple copy in case of 1 image 1114 if (nimagetot<=1) then 1115 if (use_array_all) array_img(:,:,1)=array_img_all(:,:,1) 1116 1117 else 1118 1119 ! Compute number of data 1120 rsize=size1*size2;rsize_img=nimage*rsize 1121 ABI_MALLOC(rsize_img_all,(mpi_enreg%nproc_img)) 1122 rsize_img_all(:)=rsize*nimage_all(:) 1123 1124 ! Compute shifts in buffer arrays for each proc 1125 ABI_MALLOC(rbufshft,(mpi_enreg%nproc_img)) 1126 rbufshft(1)=0 1127 do jj=2,mpi_enreg%nproc_img 1128 rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1) 1129 end do 1130 1131 ! Load buffer 1132 if (use_array_all) then 1133 ABI_MALLOC(rbuffer_all,(rsize*nimagetot)) 1134 end if 1135 if (.not.use_array_all) then 1136 ABI_MALLOC(rbuffer_all,(0)) 1137 end if 1138 if (use_array_all) then 1139 ABI_MALLOC(iimg,(mpi_enreg%nproc_img)) 1140 iimg=0 1141 do jj=1,nimagetot 1142 ! The following line supposes that images are sorted by increasing index 1143 iproc=mpi_enreg%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1 1144 ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize 1145 rbuffer_all(ibufr+1:ibufr+rsize)=reshape(array_img_all(1:size1,1:size2,jj),(/rsize/)) 1146 end do 1147 ABI_FREE(iimg) 1148 end if 1149 1150 ! Scatter all data 1151 ABI_MALLOC(rbuffer,(rsize_img)) 1152 call xmpi_scatterv(rbuffer_all,rsize_img_all,rbufshft,rbuffer,rsize_img,& 1153 & master_img,mpi_enreg%comm_img,ierr) 1154 ABI_FREE(rbuffer_all) 1155 ABI_FREE(rbufshft) 1156 ABI_FREE(rsize_img_all) 1157 1158 ! Transfered distributed buffers into array_img (master proc only) 1159 ibufr=0 1160 do jj=1,nimage 1161 array_img(1:size1,1:size2,jj)=reshape(rbuffer(ibufr+1:ibufr+rsize),(/size1,size2/)) 1162 ibufr=ibufr+rsize 1163 end do 1164 ABI_FREE(rbuffer) 1165 1166 end if ! nimagetot<=1 1167 ABI_FREE(nimage_all) 1168 1169 ! Now, is requested dispatch data inside each image 1170 if (.not.one_per_img) then 1171 call xmpi_bcast(array_img,0,mpi_enreg%comm_cell,ierr) 1172 end if 1173 1174 end if 1175 1176 end subroutine scatter_array_img