TABLE OF CONTENTS
- ABINIT/m_paw_ij
- m_paw_an/paw_ij_isendreceive_fillbuffer
- m_paw_ij/paw_ij_copy
- m_paw_ij/paw_ij_free
- m_paw_ij/paw_ij_gather
- m_paw_ij/paw_ij_init
- m_paw_ij/paw_ij_isendreceive_getbuffer
- m_paw_ij/paw_ij_nullify
- m_paw_ij/paw_ij_print
- m_paw_ij/paw_ij_redistribute
- m_paw_ij/paw_ij_reset_flags
- m_paw_ij/paw_ij_type
ABINIT/m_paw_ij [ Modules ]
NAME
m_paw_ij
FUNCTION
This module contains the definition of the paw_ij_type structured datatype, as well as related functions and methods. paw_ij_type variables contain various arrays given on (i,j) (partial waves) channels for a given atom.
COPYRIGHT
Copyright (C) 2013-2018 ABINIT group (MT, FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
24 #include "libpaw.h" 25 26 MODULE m_paw_ij 27 28 USE_DEFS 29 USE_MSG_HANDLING 30 USE_MPI_WRAPPERS 31 USE_MEMORY_PROFILING 32 33 use m_paral_atom, only : get_my_atmtab, free_my_atmtab, get_my_natom 34 use m_pawtab, only : pawtab_type 35 use m_paw_io, only : pawio_print_ij 36 37 implicit none 38 39 private
m_paw_an/paw_ij_isendreceive_fillbuffer [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_ij_isendreceive_fillbuffer
FUNCTION
Extract from paw_ij 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 paw_ij_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 npaw_ij_send= number of sent atoms paw_ij= data structure from which are extract buffer int and buffer dp
OUTPUT
buf_int= buffer of integers to be sent buf_int_size= size of buffer of integers buf_dp= buffer of double precision numbers to be sent buf_dp_size= size of buffer of double precision numbers
PARENTS
m_paw_ij
CHILDREN
SOURCE
2777 subroutine paw_ij_isendreceive_fillbuffer(paw_ij,atmtab_send,atm_indx_send,npaw_ij_send,& 2778 & buf_int,buf_int_size,buf_dp,buf_dp_size) 2779 2780 2781 !This section has been created automatically by the script Abilint (TD). 2782 !Do not modify the following lines by hand. 2783 #undef ABI_FUNC 2784 #define ABI_FUNC 'paw_ij_isendreceive_fillbuffer' 2785 !End of the abilint section 2786 2787 implicit none 2788 2789 !Arguments ------------------------------------ 2790 !scalars 2791 integer,intent(out) :: buf_int_size,buf_dp_size 2792 integer,intent(in) :: npaw_ij_send 2793 !arrays 2794 integer,intent(in) :: atmtab_send(:),atm_indx_send(:) 2795 integer,allocatable,intent(out) :: buf_int(:) 2796 real(dp),allocatable,intent(out):: buf_dp(:) 2797 type(paw_ij_type),target,intent(in) :: paw_ij(:) 2798 2799 !Local variables------------------------------- 2800 !scalars 2801 integer :: cplxdij_lmn2_size,cplxrf_lmn2_size,cplxdijrf_lmn2_size 2802 integer :: iatom_tot,ii,ij,indx_dp,indx_int,ipaw_ij_send 2803 integer :: lmn2_size,ndij,nocc,nspden,sz1,sz2,sz3 2804 character(len=500) :: msg 2805 type(Paw_ij_type),pointer :: paw_ij1 2806 !arrays 2807 2808 ! ********************************************************************* 2809 2810 !Compute sizes of buffers 2811 buf_int_size=0;buf_dp_size=0 2812 do ipaw_ij_send=1,npaw_ij_send 2813 iatom_tot=atmtab_send(ipaw_ij_send) 2814 ij = atm_indx_send(iatom_tot) 2815 paw_ij1=>paw_ij(ij) 2816 lmn2_size=paw_ij1%lmn2_size 2817 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2818 cplxrf_lmn2_size=paw_ij1%cplex_rf*lmn2_size 2819 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij1%cplex_rf 2820 ndij=paw_ij1%ndij 2821 buf_int_size=buf_int_size+24 2822 if (paw_ij1%has_dij==2) then 2823 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2824 end if 2825 if (paw_ij1%has_dij0==2) then 2826 buf_dp_size=buf_dp_size +lmn2_size 2827 end if 2828 if (paw_ij1%has_dijexxc==2) then 2829 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2830 end if 2831 if (paw_ij1%has_dijfock==2) then 2832 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2833 end if 2834 if (paw_ij1%has_dijfr==2) then 2835 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2836 end if 2837 if (paw_ij1%has_dijhartree==2) then 2838 buf_dp_size=buf_dp_size +cplxrf_lmn2_size 2839 end if 2840 if (paw_ij1%has_dijhat==2) then 2841 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2842 end if 2843 if (paw_ij1%has_dijnd==2) then 2844 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2845 end if 2846 if (paw_ij1%has_dijso==2) then 2847 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2848 end if 2849 if (paw_ij1%has_dijU==2) then 2850 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2851 end if 2852 if (paw_ij1%has_dijxc==2) then 2853 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 2854 end if 2855 if (paw_ij1%has_dijxc_hat==2) then 2856 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2857 end if 2858 if (paw_ij1%has_dijxc_val==2) then 2859 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2860 end if 2861 if (paw_ij1%has_pawu_occ>=1) then 2862 buf_int_size=buf_int_size+4 2863 buf_dp_size=buf_dp_size & 2864 & +size(paw_ij1%nocctot) & 2865 & +size(paw_ij1%noccmmp) 2866 end if 2867 if (paw_ij1%has_exexch_pot>=1) then 2868 buf_int_size=buf_int_size+3 2869 buf_dp_size=buf_dp_size +size(paw_ij1%vpawx) 2870 end if 2871 end do 2872 2873 !Fill input buffers 2874 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 2875 LIBPAW_ALLOCATE(buf_dp,(buf_dp_size)) 2876 indx_int=1;indx_dp=1 2877 do ipaw_ij_send=1,npaw_ij_send 2878 iatom_tot=atmtab_send(ipaw_ij_send) 2879 ij = atm_indx_send(iatom_tot) 2880 paw_ij1=>paw_ij(ij) 2881 nspden=paw_ij1%nspden 2882 buf_int(indx_int)=iatom_tot ;indx_int=indx_int+1 2883 buf_int(indx_int)=paw_ij1%cplex_rf ;indx_int=indx_int+1 2884 buf_int(indx_int)=paw_ij1%cplex_dij ;indx_int=indx_int+1 2885 buf_int(indx_int)=paw_ij1%itypat ;indx_int=indx_int+1 2886 buf_int(indx_int)=nspden ;indx_int=indx_int+1 2887 buf_int(indx_int)=paw_ij1%nsppol ;indx_int=indx_int+1 2888 buf_int(indx_int)=paw_ij1%lmn_size ;indx_int=indx_int+1 2889 buf_int(indx_int)=paw_ij1%lmn2_size ;indx_int=indx_int+1 2890 buf_int(indx_int)=paw_ij1%ndij ;indx_int=indx_int+1 2891 buf_int(indx_int)=paw_ij1%has_dij ;indx_int=indx_int+1 2892 buf_int(indx_int)=paw_ij1%has_dij0 ;indx_int=indx_int+1 2893 buf_int(indx_int)=paw_ij1%has_dijexxc ;indx_int=indx_int+1 2894 buf_int(indx_int)=paw_ij1%has_dijfock ;indx_int=indx_int+1 2895 buf_int(indx_int)=paw_ij1%has_dijfr ;indx_int=indx_int+1 2896 buf_int(indx_int)=paw_ij1%has_dijhartree ;indx_int=indx_int+1 2897 buf_int(indx_int)=paw_ij1%has_dijhat ;indx_int=indx_int+1 2898 buf_int(indx_int)=paw_ij1%has_dijnd ;indx_int=indx_int+1 2899 buf_int(indx_int)=paw_ij1%has_dijso ;indx_int=indx_int+1 2900 buf_int(indx_int)=paw_ij1%has_dijU ;indx_int=indx_int+1 2901 buf_int(indx_int)=paw_ij1%has_dijxc ;indx_int=indx_int+1 2902 buf_int(indx_int)=paw_ij1%has_dijxc_hat ;indx_int=indx_int+1 2903 buf_int(indx_int)=paw_ij1%has_dijxc_val ;indx_int=indx_int+1 2904 buf_int(indx_int)=paw_ij1%has_exexch_pot ;indx_int=indx_int+1 2905 buf_int(indx_int)=paw_ij1%has_pawu_occ ;indx_int=indx_int+1 2906 lmn2_size=paw_ij1%lmn2_size 2907 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2908 cplxrf_lmn2_size=paw_ij1%cplex_rf*lmn2_size 2909 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij1%cplex_rf 2910 ndij=paw_ij1%ndij 2911 if (paw_ij1%has_dij==2) then 2912 ii=cplxdijrf_lmn2_size*ndij 2913 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dij,(/ii/)) 2914 indx_dp=indx_dp+ii 2915 end if 2916 if (paw_ij1%has_dij0==2) then 2917 ii=lmn2_size 2918 buf_dp(indx_dp:indx_dp+lmn2_size-1)=paw_ij1%dij0(:) 2919 indx_dp=indx_dp+lmn2_size 2920 end if 2921 if (paw_ij1%has_dijexxc==2) then 2922 ii=cplxdij_lmn2_size*ndij 2923 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijexxc,(/ii/)) 2924 indx_dp=indx_dp+ii 2925 end if 2926 if (paw_ij1%has_dijfock==2) then 2927 ii=cplxdij_lmn2_size*ndij 2928 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijfock,(/ii/)) 2929 indx_dp=indx_dp+ii 2930 end if 2931 if (paw_ij1%has_dijfr==2) then 2932 ii=cplxdijrf_lmn2_size*ndij 2933 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijfr,(/ii/)) 2934 indx_dp=indx_dp+ii 2935 end if 2936 if (paw_ij1%has_dijhartree==2) then 2937 ii=cplxrf_lmn2_size 2938 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij1%dijhartree(:) 2939 indx_dp=indx_dp+ii 2940 end if 2941 if (paw_ij1%has_dijhat==2) then 2942 ii=cplxdijrf_lmn2_size*ndij 2943 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijhat,(/ii/)) 2944 indx_dp=indx_dp+ii 2945 end if 2946 if (paw_ij1%has_dijnd==2) then 2947 ii=cplxdij_lmn2_size*ndij 2948 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijnd,(/ii/)) 2949 indx_dp=indx_dp+ii 2950 end if 2951 if (paw_ij1%has_dijso==2) then 2952 ii=cplxdijrf_lmn2_size*ndij 2953 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijso,(/ii/)) 2954 indx_dp=indx_dp+ii 2955 end if 2956 if (paw_ij1%has_dijU==2) then 2957 ii=cplxdijrf_lmn2_size*ndij 2958 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijU,(/ii/)) 2959 indx_dp=indx_dp+ii 2960 end if 2961 if (paw_ij1%has_dijxc==2) then 2962 ii=cplxdijrf_lmn2_size*ndij 2963 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc,(/ii/)) 2964 indx_dp=indx_dp+ii 2965 end if 2966 if (paw_ij1%has_dijxc_hat==2) then 2967 ii=cplxdij_lmn2_size*ndij 2968 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc_hat,(/ii/)) 2969 indx_dp=indx_dp+ii 2970 end if 2971 if (paw_ij1%has_dijxc_val==2) then 2972 ii=cplxdij_lmn2_size*ndij 2973 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc_val,(/ii/)) 2974 indx_dp=indx_dp+ii 2975 end if 2976 if (paw_ij1%has_pawu_occ>=1)then 2977 buf_int(indx_int)=size(paw_ij1%noccmmp,1) ;indx_int=indx_int+1 2978 buf_int(indx_int)=size(paw_ij1%noccmmp,2) ;indx_int=indx_int+1 2979 buf_int(indx_int)=size(paw_ij1%noccmmp,3) ;indx_int=indx_int+1 2980 buf_int(indx_int)=size(paw_ij1%noccmmp,4) ;indx_int=indx_int+1 2981 nocc=paw_ij1%ndij 2982 buf_dp(indx_dp:indx_dp+nocc-1)=paw_ij1%nocctot(1:nocc) 2983 indx_dp=indx_dp+nocc 2984 nocc=size(paw_ij1%noccmmp) 2985 buf_dp(indx_dp:indx_dp+nocc-1)=reshape(paw_ij1%noccmmp,(/nocc/)) 2986 indx_dp=indx_dp+nocc 2987 end if 2988 if (paw_ij1%has_exexch_pot>=1) then 2989 sz1=size(paw_ij1%vpawx,1);sz2=size(paw_ij1%vpawx,2) 2990 sz3=size(paw_ij1%vpawx,3) 2991 buf_int(indx_int)=sz1; indx_int=indx_int+1 2992 buf_int(indx_int)=sz2; indx_int=indx_int+1 2993 buf_int(indx_int)=sz3; indx_int=indx_int+1 2994 ii=sz1*sz2*sz3 2995 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%vpawx,(/(ii)/)) 2996 indx_dp=indx_dp+ii 2997 end if 2998 end do 2999 indx_int=indx_int-1;indx_dp=indx_dp-1 3000 if ((indx_int/=buf_int_size).or.(indx_dp/=buf_dp_size)) then 3001 write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size 3002 MSG_BUG(msg) 3003 end if 3004 3005 end subroutine paw_ij_isendreceive_fillbuffer
m_paw_ij/paw_ij_copy [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_copy
FUNCTION
Copy one paw_ij datastructure into another Can take into accound changes of dimensions Can copy a shared paw_ij 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 paw_ij_in(:)<type(paw_ij_type)>= input paw_ij datastructure
SIDE EFFECTS
paw_ij_cpy(:)<type(paw_ij_type)>= output paw_ij datastructure
NOTES
paw_ij_cpy must have been allocated in the calling function.
PARENTS
m_paw_ij
CHILDREN
SOURCE
775 subroutine paw_ij_copy(paw_ij_in,paw_ij_cpy, & 776 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 777 778 779 !This section has been created automatically by the script Abilint (TD). 780 !Do not modify the following lines by hand. 781 #undef ABI_FUNC 782 #define ABI_FUNC 'paw_ij_copy' 783 !End of the abilint section 784 785 implicit none 786 787 !Arguments ------------------------------------ 788 !scalars 789 integer,optional,intent(in) :: comm_atom 790 !arrays 791 integer,optional,target,intent(in) :: mpi_atmtab(:) 792 type(paw_ij_type),intent(in) :: paw_ij_in(:) 793 type(paw_ij_type),intent(inout),target :: paw_ij_cpy(:) 794 795 !Local variables------------------------------- 796 !scalars 797 integer :: ij,ij1,my_comm_atom,my_paw_ij,npaw_ij_in,npaw_ij_max,npaw_ij_out,paral_case,sz1,sz2,sz3,sz4 798 logical :: my_atmtab_allocated,paral_atom 799 character(len=500) :: msg 800 !arrays 801 integer,pointer :: my_atmtab(:) 802 type(paw_ij_type),pointer :: paw_ij_out(:) 803 804 ! ************************************************************************* 805 806 !@Paw_ij_type 807 808 !Retrieve sizes 809 npaw_ij_in=size(paw_ij_in);npaw_ij_out=size(paw_ij_cpy) 810 811 !Set up parallelism over atoms 812 paral_atom=(present(comm_atom));if (paral_atom) paral_atom=(xmpi_comm_size(comm_atom)>1) 813 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 814 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 815 my_atmtab_allocated=.false. 816 817 !Determine in which case we are (parallelism, ...) 818 !No parallelism: a single copy operation 819 paral_case=0;npaw_ij_max=npaw_ij_in 820 paw_ij_out => paw_ij_cpy 821 if (paral_atom) then 822 if (npaw_ij_out<npaw_ij_in) then ! Parallelism: the copy operation is a scatter 823 call get_my_natom(my_comm_atom,my_paw_ij,npaw_ij_in) 824 if (my_paw_ij==npaw_ij_out) then 825 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,npaw_ij_in) 826 paral_case=1;npaw_ij_max=npaw_ij_out 827 paw_ij_out => paw_ij_cpy 828 else 829 msg=' npaw_ij_out should be equal to my_paw_ij !' 830 MSG_BUG(msg) 831 end if 832 else ! Parallelism: the copy operation is a gather 833 call get_my_natom(my_comm_atom,my_paw_ij,npaw_ij_out) 834 if (my_paw_ij==npaw_ij_in) then 835 paral_case=2;npaw_ij_max=npaw_ij_in 836 else 837 msg=' npaw_ij_in should be equal to my_paw_ij !' 838 MSG_BUG(msg) 839 end if 840 end if 841 end if 842 843 !First case: a simple copy or a scatter 844 if (npaw_ij_max>0.and.((paral_case==0).or.(paral_case==1))) then 845 call paw_ij_nullify(paw_ij_out) 846 847 ! Loop on paw_ij components 848 do ij1=1,npaw_ij_max 849 ij=ij1; if (paral_case==1) ij=my_atmtab(ij1) 850 851 paw_ij_out(ij1)%cplex_rf=paw_ij_in(ij)%cplex_rf 852 paw_ij_out(ij1)%cplex_dij=paw_ij_in(ij)%cplex_dij 853 paw_ij_out(ij1)%has_dij=paw_ij_in(ij)%has_dij 854 paw_ij_out(ij1)%has_dij0=paw_ij_in(ij)%has_dij0 855 paw_ij_out(ij1)%has_dijexxc=paw_ij_in(ij)%has_dijexxc 856 paw_ij_out(ij1)%has_dijfock=paw_ij_in(ij)%has_dijfock 857 paw_ij_out(ij1)%has_dijfr=paw_ij_in(ij)%has_dijfr 858 paw_ij_out(ij1)%has_dijhartree=paw_ij_in(ij)%has_dijhartree 859 paw_ij_out(ij1)%has_dijhat=paw_ij_in(ij)%has_dijhat 860 paw_ij_out(ij1)%has_dijnd=paw_ij_in(ij)%has_dijnd 861 paw_ij_out(ij1)%has_dijso=paw_ij_in(ij)%has_dijso 862 paw_ij_out(ij1)%has_dijU=paw_ij_in(ij)%has_dijU 863 paw_ij_out(ij1)%has_dijxc=paw_ij_in(ij)%has_dijxc 864 paw_ij_out(ij1)%has_dijxc_hat=paw_ij_in(ij)%has_dijxc_hat 865 paw_ij_out(ij1)%has_dijxc_val=paw_ij_in(ij)%has_dijxc_val 866 paw_ij_out(ij1)%has_exexch_pot=paw_ij_in(ij)%has_exexch_pot 867 paw_ij_out(ij1)%has_pawu_occ=paw_ij_in(ij)%has_pawu_occ 868 paw_ij_out(ij1)%itypat=paw_ij_in(ij)%itypat 869 paw_ij_out(ij1)%lmn_size=paw_ij_in(ij)%lmn_size 870 paw_ij_out(ij1)%lmn2_size=paw_ij_in(ij)%lmn2_size 871 paw_ij_out(ij1)%ndij=paw_ij_in(ij)%ndij 872 paw_ij_out(ij1)%nspden=paw_ij_in(ij)%nspden 873 paw_ij_out(ij1)%nsppol=paw_ij_in(ij)%nsppol 874 if (paw_ij_in(ij)%has_dij>=1) then 875 sz1=size(paw_ij_in(ij)%dij,1);sz2=size(paw_ij_in(ij)%dij,2) 876 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dij,(sz1,sz2)) 877 if (paw_ij_in(ij)%has_dij==2) & 878 & paw_ij_out(ij1)%dij(:,:)=paw_ij_in(ij)%dij(:,:) 879 end if 880 if (paw_ij_in(ij)%has_dij0>=1) then 881 sz1=size(paw_ij_in(ij)%dij0,1) 882 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dij0,(sz1)) 883 if (paw_ij_in(ij)%has_dij0 ==2) & 884 & paw_ij_out(ij1)%dij0(:)=paw_ij_in(ij)%dij0(:) 885 end if 886 if (paw_ij_in(ij)%has_dijexxc>=1) then 887 sz1=size(paw_ij_in(ij)%dijexxc,1);sz2=size(paw_ij_in(ij)%dijexxc,2) 888 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijexxc,(sz1,sz2)) 889 if (paw_ij_in(ij)%has_dijexxc==2) & 890 & paw_ij_out(ij1)%dijexxc(:,:)=paw_ij_in(ij)%dijexxc(:,:) 891 end if 892 if (paw_ij_in(ij)%has_dijfock>=1) then 893 sz1=size(paw_ij_in(ij)%dijfock,1);sz2=size(paw_ij_in(ij)%dijfock,2) 894 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijfock,(sz1,sz2)) 895 if (paw_ij_in(ij)%has_dijfock==2) & 896 & paw_ij_out(ij1)%dijfock(:,:)=paw_ij_in(ij)%dijfock(:,:) 897 end if 898 if (paw_ij_in(ij)%has_dijfr>=1) then 899 sz1=size(paw_ij_in(ij)%dijfr,1);sz2=size(paw_ij_in(ij)%dijfr,2) 900 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijfr,(sz1,sz2)) 901 if (paw_ij_in(ij)%has_dijfr==2) & 902 & paw_ij_out(ij1)%dijfr(:,:)=paw_ij_in(ij)%dijfr(:,:) 903 end if 904 if (paw_ij_in(ij)%has_dijhartree>=1) then 905 sz1=size(paw_ij_in(ij)%dijhartree,1) 906 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijhartree,(sz1)) 907 if (paw_ij_in(ij)%has_dijhartree==2) & 908 & paw_ij_out(ij1)%dijhartree(:)=paw_ij_in(ij)%dijhartree(:) 909 end if 910 if (paw_ij_in(ij)%has_dijhat>=1) then 911 sz1=size(paw_ij_in(ij)%dijhat,1);sz2=size(paw_ij_in(ij)%dijhat,2) 912 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijhat,(sz1,sz2)) 913 if (paw_ij_in(ij)%has_dijhat==2) & 914 & paw_ij_out(ij1)%dijhat(:,:)=paw_ij_in(ij)%dijhat(:,:) 915 end if 916 if (paw_ij_in(ij)%has_dijnd>=1) then 917 sz1=size(paw_ij_in(ij)%dijnd,1);sz2=size(paw_ij_in(ij)%dijnd,2) 918 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijnd,(sz1,sz2)) 919 if (paw_ij_in(ij)%has_dijnd==2) & 920 & paw_ij_out(ij1)%dijnd(:,:)=paw_ij_in(ij)%dijnd(:,:) 921 end if 922 if (paw_ij_in(ij)%has_dijU>=1) then 923 sz1=size(paw_ij_in(ij)%dijU,1);sz2=size(paw_ij_in(ij)%dijU,2) 924 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijU,(sz1,sz2)) 925 if (paw_ij_in(ij)%has_dijU==2) & 926 & paw_ij_out(ij1)%dijU(:,:)=paw_ij_in(ij)%dijU(:,:) 927 end if 928 if (paw_ij_in(ij)%has_dijso>=1) then 929 sz1=size(paw_ij_in(ij)%dijso,1);sz2=size(paw_ij_in(ij)%dijso,2) 930 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijso,(sz1,sz2)) 931 if (paw_ij_in(ij)%has_dijso==2) & 932 & paw_ij_out(ij1)%dijso(:,:)=paw_ij_in(ij)%dijso(:,:) 933 end if 934 if (paw_ij_in(ij)%has_dijxc>=1) then 935 sz1=size(paw_ij_in(ij)%dijxc,1);sz2=size(paw_ij_in(ij)%dijxc,2) 936 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc,(sz1,sz2)) 937 if (paw_ij_in(ij)%has_dijxc==2) & 938 & paw_ij_out(ij1)%dijxc(:,:)=paw_ij_in(ij)%dijxc(:,:) 939 end if 940 if (paw_ij_in(ij)%has_dijxc_hat>=1) then 941 sz1=size(paw_ij_in(ij)%dijxc_hat,1);sz2=size(paw_ij_in(ij)%dijxc_hat,2) 942 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc_hat,(sz1,sz2)) 943 if (paw_ij_in(ij)%has_dijxc_hat==2) & 944 & paw_ij_out(ij1)%dijxc_hat(:,:)=paw_ij_in(ij)%dijxc_hat(:,:) 945 end if 946 if (paw_ij_in(ij)%has_dijxc_val>=1) then 947 sz1=size(paw_ij_in(ij)%dijxc_val,1);sz2=size(paw_ij_in(ij)%dijxc_val,2) 948 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc_val,(sz1,sz2)) 949 if (paw_ij_in(ij)%has_dijxc_val==2) & 950 & paw_ij_out(ij1)%dijxc_val(:,:)=paw_ij_in(ij)%dijxc_val(:,:) 951 end if 952 if (paw_ij_in(ij)%has_pawu_occ>=1) then 953 sz1=size(paw_ij_in(ij)%noccmmp,1);sz2=size(paw_ij_in(ij)%noccmmp,2) 954 sz3=size(paw_ij_in(ij)%noccmmp,3);sz4=size(paw_ij_in(ij)%noccmmp,4) 955 LIBPAW_ALLOCATE(paw_ij_out(ij1)%noccmmp,(sz1,sz2,sz3,sz4)) 956 sz1=size(paw_ij_in(ij)%nocctot,1) 957 LIBPAW_ALLOCATE(paw_ij_out(ij1)%nocctot,(sz1)) 958 if (paw_ij_in(ij)%has_pawu_occ==2) then 959 paw_ij_out(ij1)%noccmmp(:,:,:,:)=paw_ij_in(ij)%noccmmp(:,:,:,:) 960 paw_ij_out(ij1)%nocctot(:)=paw_ij_in(ij)%nocctot(:) 961 end if 962 end if 963 if (paw_ij_in(ij)%has_exexch_pot >= 1) then 964 sz1=size(paw_ij_in(ij)%vpawx,1);sz2=size(paw_ij_in(ij)%vpawx,2) 965 sz3=size(paw_ij_in(ij)%vpawx,3) 966 LIBPAW_ALLOCATE(paw_ij_out(ij1)%vpawx,(sz1,sz2,sz3)) 967 if (paw_ij_in(ij)%has_exexch_pot==2) & 968 & paw_ij_out(ij1)%vpawx(:,:,:)=paw_ij_in(ij)%vpawx(:,:,:) 969 end if 970 971 end do 972 end if 973 974 !Second case: a gather 975 if (paral_case==2) then 976 call paw_ij_gather(paw_ij_in,paw_ij_cpy,-1,my_comm_atom) 977 end if 978 979 !Destroy atom table used for parallelism 980 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 981 982 end subroutine paw_ij_copy
m_paw_ij/paw_ij_free [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_free
FUNCTION
Deallocate pointers and nullify flags in a paw_ij structure
SIDE EFFECTS
paw_ij(:)<type(paw_ij_type)>=paw arrays given on (i,j) channels
PARENTS
bethe_salpeter,d2frnl,dfpt_nstpaw,dfpt_rhofermi,dfpt_scfcv,ldau_self m_energy,m_paral_pert,m_paw_ij,pawprt,respfn,scfcv,screening,sigma
CHILDREN
SOURCE
578 subroutine paw_ij_free(Paw_ij) 579 580 581 !This section has been created automatically by the script Abilint (TD). 582 !Do not modify the following lines by hand. 583 #undef ABI_FUNC 584 #define ABI_FUNC 'paw_ij_free' 585 !End of the abilint section 586 587 implicit none 588 589 !Arguments ------------------------------------ 590 !arrays 591 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 592 593 !Local variables------------------------------- 594 integer :: iat,natom 595 596 ! ************************************************************************* 597 598 !@Paw_ij_type 599 600 natom=SIZE(Paw_ij);if (natom==0) return 601 602 do iat=1,natom 603 if (allocated(Paw_ij(iat)%dij )) then 604 LIBPAW_DEALLOCATE(Paw_ij(iat)%dij) 605 end if 606 if (allocated(Paw_ij(iat)%dij0 )) then 607 LIBPAW_DEALLOCATE(Paw_ij(iat)%dij0) 608 end if 609 if (allocated(Paw_ij(iat)%dijexxc )) then 610 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijexxc) 611 end if 612 if (allocated(Paw_ij(iat)%dijfock )) then 613 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijfock) 614 end if 615 if (allocated(Paw_ij(iat)%dijfr )) then 616 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijfr) 617 end if 618 if (allocated(Paw_ij(iat)%dijhartree)) then 619 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijhartree) 620 end if 621 if (allocated(Paw_ij(iat)%dijhat )) then 622 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijhat) 623 end if 624 if (allocated(Paw_ij(iat)%dijnd )) then 625 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijnd) 626 end if 627 if (allocated(Paw_ij(iat)%dijU )) then 628 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijU) 629 end if 630 if (allocated(Paw_ij(iat)%dijso )) then 631 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijso) 632 end if 633 if (allocated(Paw_ij(iat)%dijxc )) then 634 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc) 635 end if 636 if (allocated(Paw_ij(iat)%dijxc_hat )) then 637 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc_hat) 638 end if 639 if (allocated(Paw_ij(iat)%dijxc_val )) then 640 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc_val) 641 end if 642 if (allocated(Paw_ij(iat)%noccmmp )) then 643 LIBPAW_DEALLOCATE(Paw_ij(iat)%noccmmp) 644 end if 645 if (allocated(Paw_ij(iat)%nocctot )) then 646 LIBPAW_DEALLOCATE(Paw_ij(iat)%nocctot) 647 end if 648 if (allocated(Paw_ij(iat)%vpawx )) then 649 LIBPAW_DEALLOCATE(Paw_ij(iat)%vpawx) 650 end if 651 652 ! === Reset all has_* flags === 653 Paw_ij(iat)%has_dij =0 654 Paw_ij(iat)%has_dij0 =0 655 Paw_ij(iat)%has_dijexxc =0 656 Paw_ij(iat)%has_dijfock =0 657 Paw_ij(iat)%has_dijfr =0 658 Paw_ij(iat)%has_dijhartree=0 659 Paw_ij(iat)%has_dijhat =0 660 Paw_ij(iat)%has_dijnd =0 661 Paw_ij(iat)%has_dijso =0 662 Paw_ij(iat)%has_dijU =0 663 Paw_ij(iat)%has_dijxc =0 664 Paw_ij(iat)%has_dijxc_hat =0 665 Paw_ij(iat)%has_dijxc_val =0 666 Paw_ij(iat)%has_exexch_pot=0 667 Paw_ij(iat)%has_pawu_occ =0 668 end do 669 670 !call paw_ij_nullify(Paw_ij) 671 672 end subroutine paw_ij_free
m_paw_ij/paw_ij_gather [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_gather
FUNCTION
(All)Gather paw_ij datastructures
INPUTS
master=master receiving data ; if -1 do a ALLGATHER comm_atom= communicator paw_ij_in(:)<type(paw_ij_type)>= input paw_ij datastructures on every process
OUTPUT
paw_ij_gathered(:)<type(paw_ij_type)>= output paw_oij datastructure
PARENTS
m_paw_ij,outkss,pawprt
CHILDREN
SOURCE
1404 subroutine paw_ij_gather(paw_ij_in,paw_ij_gathered,master,comm_atom) 1405 1406 1407 !This section has been created automatically by the script Abilint (TD). 1408 !Do not modify the following lines by hand. 1409 #undef ABI_FUNC 1410 #define ABI_FUNC 'paw_ij_gather' 1411 !End of the abilint section 1412 1413 implicit none 1414 1415 !Arguments ------------------------------------ 1416 !scalars 1417 integer,intent(in) :: master,comm_atom 1418 !arrays 1419 type(paw_ij_type),intent(in) :: paw_ij_in(:) 1420 type(paw_ij_type),intent(inout) :: paw_ij_gathered(:) 1421 1422 !Local variables------------------------------- 1423 !scalars 1424 integer :: buf_dp_size,buf_dp_size_all,buf_int_size,buf_int_size_all 1425 integer :: cplxrf_lmn2_size,cplxdij_lmn2_size,cplxdijrf_lmn2_size 1426 integer :: iat,ii,ierr,ij,indx_dp,indx_int,lmn2_size 1427 integer :: me_atom,ndij,nocc,nocc1,nocc2,nocc3,nocc4,npaw_ij_in,npaw_ij_in_sum 1428 integer :: npaw_ij_out,nproc_atom,nspden,sz1,sz2,sz3,sz4 1429 logical :: my_atmtab_allocated,paral_atom 1430 character(len=500) :: msg 1431 !arrays 1432 integer :: bufsz(2) 1433 integer,allocatable :: buf_int(:),buf_int_all(:) 1434 integer,allocatable :: count_dp(:),count_int(:),count_tot(:),displ_dp(:),displ_int(:) 1435 integer,pointer :: my_atmtab(:) 1436 real(dp),allocatable :: buf_dp(:),buf_dp_all(:) 1437 1438 ! ************************************************************************* 1439 1440 !@Paw_ij_type 1441 1442 npaw_ij_in=size(paw_ij_in);npaw_ij_out=size(paw_ij_gathered) 1443 1444 nproc_atom=xmpi_comm_size(comm_atom) 1445 me_atom=xmpi_comm_rank(comm_atom) 1446 1447 !Special case 1 process 1448 if (nproc_atom==1) then 1449 if (master==-1.or.me_atom==master) then 1450 call paw_ij_free(paw_ij_gathered) 1451 call paw_ij_nullify(paw_ij_gathered) 1452 do iat=1,npaw_ij_in 1453 paw_ij_gathered(iat)%cplex_rf =paw_ij_in(iat)%cplex_rf 1454 paw_ij_gathered(iat)%cplex_dij =paw_ij_in(iat)%cplex_dij 1455 Paw_ij_gathered(iat)%has_dij =paw_ij_in(iat)%has_dij 1456 Paw_ij_gathered(iat)%has_dij0 =paw_ij_in(iat)%has_dij0 1457 Paw_ij_gathered(iat)%has_dijexxc =paw_ij_in(iat)%has_dijexxc 1458 Paw_ij_gathered(iat)%has_dijfock =paw_ij_in(iat)%has_dijfock 1459 Paw_ij_gathered(iat)%has_dijfr =paw_ij_in(iat)%has_dijfr 1460 Paw_ij_gathered(iat)%has_dijhartree=paw_ij_in(iat)%has_dijhartree 1461 Paw_ij_gathered(iat)%has_dijhat =paw_ij_in(iat)%has_dijhat 1462 Paw_ij_gathered(iat)%has_dijnd =paw_ij_in(iat)%has_dijnd 1463 Paw_ij_gathered(iat)%has_dijso =paw_ij_in(iat)%has_dijso 1464 Paw_ij_gathered(iat)%has_dijU =paw_ij_in(iat)%has_dijU 1465 Paw_ij_gathered(iat)%has_dijxc =paw_ij_in(iat)%has_dijxc 1466 Paw_ij_gathered(iat)%has_dijxc_hat =paw_ij_in(iat)%has_dijxc_hat 1467 Paw_ij_gathered(iat)%has_dijxc_val =paw_ij_in(iat)%has_dijxc_val 1468 Paw_ij_gathered(iat)%has_exexch_pot=paw_ij_in(iat)%has_exexch_pot 1469 Paw_ij_gathered(iat)%has_pawu_occ =paw_ij_in(iat)%has_pawu_occ 1470 paw_ij_gathered(iat)%itypat =paw_ij_in(iat)%itypat 1471 paw_ij_gathered(iat)%lmn_size =paw_ij_in(iat)%lmn_size 1472 paw_ij_gathered(iat)%lmn2_size =paw_ij_in(iat)%lmn2_size 1473 paw_ij_gathered(iat)%ndij =paw_ij_in(iat)%ndij 1474 paw_ij_gathered(iat)%nspden =paw_ij_in(iat)%nspden 1475 paw_ij_gathered(iat)%nsppol =paw_ij_in(iat)%nsppol 1476 lmn2_size=paw_ij_gathered(iat)%lmn2_size 1477 cplxdij_lmn2_size=paw_ij_gathered(iat)%cplex_dij*lmn2_size 1478 cplxrf_lmn2_size=paw_ij_gathered(iat)%cplex_rf*lmn2_size 1479 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij_gathered(iat)%cplex_rf 1480 ndij=paw_ij_gathered(iat)%ndij 1481 nspden=paw_ij_gathered(iat)%nspden 1482 if (paw_ij_gathered(iat)%has_dij>=1) then 1483 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij,(cplxdijrf_lmn2_size,ndij)) 1484 if (paw_ij_in(iat)%has_dij==2) then 1485 paw_ij_gathered(iat)%dij=paw_ij_in(iat)%dij 1486 end if 1487 end if 1488 if (paw_ij_gathered(iat)%has_dij0 >=1) then 1489 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij0,(lmn2_size)) 1490 if (paw_ij_in(iat)%has_dij0==2) then 1491 paw_ij_gathered(iat)%dij0=paw_ij_in(iat)%dij0 1492 end if 1493 end if 1494 if (paw_ij_gathered(iat)%has_dijexxc >=1) then 1495 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijexxc,(cplxdij_lmn2_size,ndij)) 1496 if (paw_ij_in(iat)%has_dijexxc==2) then 1497 paw_ij_gathered(iat)%dijexxc(:,:)=paw_ij_in(iat)%dijexxc(:,:) 1498 end if 1499 end if 1500 if (paw_ij_gathered(iat)%has_dijfock >=1) then 1501 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfock,(cplxdij_lmn2_size,ndij)) 1502 if (paw_ij_in(iat)%has_dijfock==2) then 1503 paw_ij_gathered(iat)%dijfock(:,:)=paw_ij_in(iat)%dijfock(:,:) 1504 end if 1505 end if 1506 if (paw_ij_gathered(iat)%has_dijfr >=1) then 1507 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfr,(cplxdijrf_lmn2_size,ndij)) 1508 if (paw_ij_in(iat)%has_dijfr==2) then 1509 paw_ij_gathered(iat)%dijfr(:,:)=paw_ij_in(iat)%dijfr(:,:) 1510 end if 1511 end if 1512 if (paw_ij_gathered(iat)%has_dijhartree >=1) then 1513 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhartree,(cplxrf_lmn2_size)) 1514 if (paw_ij_in(iat)%has_dijhartree==2) then 1515 paw_ij_gathered(iat)%dijhartree(:)=paw_ij_in(iat)%dijhartree(:) 1516 end if 1517 end if 1518 if (paw_ij_gathered(iat)%has_dijhat >=1) then 1519 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhat,(cplxdijrf_lmn2_size,ndij)) 1520 if (paw_ij_in(iat)%has_dijhat==2) then 1521 paw_ij_gathered(iat)%dijhat(:,:)=paw_ij_in(iat)%dijhat(:,:) 1522 end if 1523 end if 1524 if (paw_ij_gathered(iat)%has_dijnd >=1) then 1525 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijnd,(cplxdij_lmn2_size,ndij)) 1526 if (paw_ij_in(iat)%has_dijnd==2) then 1527 paw_ij_gathered(iat)%dijnd(:,:)=paw_ij_in(iat)%dijnd(:,:) 1528 end if 1529 end if 1530 if (paw_ij_gathered(iat)%has_dijU >=1) then 1531 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijU,(cplxdijrf_lmn2_size,ndij)) 1532 if (paw_ij_in(iat)%has_dijU==2) then 1533 paw_ij_gathered(iat)%dijU(:,:)=paw_ij_in(iat)%dijU(:,:) 1534 end if 1535 end if 1536 if (paw_ij_gathered(iat)%has_dijso >=1) then 1537 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijso,(cplxdijrf_lmn2_size,ndij)) 1538 if (paw_ij_in(iat)%has_dijso==2) then 1539 paw_ij_gathered(iat)%dijso(:,:)=paw_ij_in(iat)%dijso(:,:) 1540 end if 1541 end if 1542 if (paw_ij_gathered(iat)%has_dijxc >=1) then 1543 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc,(cplxdijrf_lmn2_size,ndij)) 1544 if (paw_ij_in(iat)%has_dijxc==2) then 1545 paw_ij_gathered(iat)%dijxc(:,:)=paw_ij_in(iat)%dijxc(:,:) 1546 end if 1547 end if 1548 if (paw_ij_gathered(iat)%has_dijxc_hat >=1) then 1549 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_hat,(cplxdij_lmn2_size,ndij)) 1550 if (paw_ij_in(iat)%has_dijxc_hat==2) then 1551 paw_ij_gathered(iat)%dijxc_hat(:,:)=paw_ij_in(iat)%dijxc_hat(:,:) 1552 end if 1553 end if 1554 if (paw_ij_gathered(iat)%has_dijxc_val >=1) then 1555 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_val,(cplxdij_lmn2_size,ndij)) 1556 if (paw_ij_in(iat)%has_dijxc_val==2) then 1557 paw_ij_gathered(iat)%dijxc_val(:,:)=paw_ij_in(iat)%dijxc_val(:,:) 1558 end if 1559 end if 1560 if (paw_ij_gathered(iat)%has_pawu_occ >=1) then 1561 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%nocctot,(ndij)) 1562 if (paw_ij_gathered(iat)%has_pawu_occ==2) then 1563 paw_ij_gathered(iat)%nocctot(:)=paw_ij_in(iat)%nocctot(:) 1564 if (allocated(paw_ij_in(iat)%noccmmp)) then 1565 sz1=size(paw_ij_in(iat)%noccmmp,1);sz2=size(paw_ij_in(iat)%noccmmp,2) 1566 sz3=size(paw_ij_in(iat)%noccmmp,3);sz4=size(paw_ij_in(iat)%noccmmp,4) 1567 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%noccmmp,(sz1,sz2,sz3,sz4)) 1568 paw_ij_gathered(iat)%noccmmp(:,:,:,:)=paw_ij_in(iat)%noccmmp(:,:,:,:) 1569 end if 1570 end if 1571 end if 1572 if (paw_ij_in(iat)%has_exexch_pot >=1) then 1573 sz1=size(paw_ij_in(iat)%vpawx,1);sz2=size(paw_ij_in(iat)%vpawx,2) 1574 sz3=size(paw_ij_in(iat)%vpawx,3) 1575 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%vpawx,(sz1,sz2,sz3)) 1576 if (paw_ij_in(iat)%has_exexch_pot==2) then 1577 paw_ij_gathered(iat)%vpawx(:,:,:)=paw_ij_in(iat)%vpawx(:,:,:) 1578 end if 1579 end if 1580 end do ! iat 1581 end if 1582 return 1583 end if !nproc_atom =1 1584 1585 !Test on sizes 1586 npaw_ij_in_sum=npaw_ij_in 1587 1588 call xmpi_sum(npaw_ij_in_sum,comm_atom,ierr) 1589 if (master==-1) then 1590 if (npaw_ij_out/=npaw_ij_in_sum) then 1591 msg='Wrong sizes sum[npaw_ij_ij]/=npaw_ij_out !' 1592 MSG_BUG(msg) 1593 end if 1594 else 1595 if (me_atom==master.and.npaw_ij_out/=npaw_ij_in_sum) then 1596 msg='(2) paw_ij_gathered wrongly allocated !' 1597 MSG_BUG(msg) 1598 end if 1599 end if 1600 1601 !Retrieve table of atoms 1602 paral_atom=.true.;nullify(my_atmtab) 1603 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,npaw_ij_in_sum) 1604 1605 !Compute sizes of buffers 1606 buf_int_size=0;buf_dp_size=0 1607 do ij=1,npaw_ij_in 1608 lmn2_size=paw_ij_in(ij)%lmn2_size 1609 cplxdij_lmn2_size=paw_ij_in(ij)%cplex_dij*lmn2_size 1610 cplxrf_lmn2_size=paw_ij_in(ij)%cplex_rf*lmn2_size 1611 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij_in(ij)%cplex_rf 1612 ndij=paw_ij_in(ij)%ndij 1613 buf_int_size=buf_int_size+24 1614 if (paw_ij_in(ij)%has_dij==2) then 1615 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1616 end if 1617 if (paw_ij_in(ij)%has_dij0==2) then 1618 buf_dp_size=buf_dp_size +lmn2_size 1619 end if 1620 if (paw_ij_in(ij)%has_dijexxc==2) then 1621 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1622 end if 1623 if (paw_ij_in(ij)%has_dijfock==2) then 1624 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1625 end if 1626 if (paw_ij_in(ij)%has_dijfr==2) then 1627 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1628 end if 1629 if (paw_ij_in(ij)%has_dijhartree==2) then 1630 buf_dp_size=buf_dp_size +cplxrf_lmn2_size 1631 end if 1632 if (paw_ij_in(ij)%has_dijhat==2) then 1633 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1634 end if 1635 if (paw_ij_in(ij)%has_dijnd==2) then 1636 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1637 end if 1638 if (paw_ij_in(ij)%has_dijso==2) then 1639 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1640 end if 1641 if (paw_ij_in(ij)%has_dijU==2) then 1642 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1643 end if 1644 if (paw_ij_in(ij)%has_dijxc==2) then 1645 buf_dp_size=buf_dp_size +cplxdijrf_lmn2_size*ndij 1646 end if 1647 if (paw_ij_in(ij)%has_dijxc_hat==2) then 1648 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1649 end if 1650 if (paw_ij_in(ij)%has_dijxc_val==2) then 1651 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1652 end if 1653 if (paw_ij_in(ij)%has_pawu_occ>=1) then 1654 buf_int_size=buf_int_size+4 1655 buf_dp_size=buf_dp_size & 1656 & +size(paw_ij_in(ij)%nocctot) & 1657 & +size(paw_ij_in(ij)%noccmmp) 1658 end if 1659 if (paw_ij_in(ij)%has_exexch_pot>=1) then 1660 buf_int_size=buf_int_size+3 1661 buf_dp_size=buf_dp_size +size(paw_ij_in(ij)%vpawx) 1662 end if 1663 end do 1664 1665 !Fill input buffers 1666 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 1667 LIBPAW_ALLOCATE(buf_dp ,(buf_dp_size)) 1668 indx_int=1;indx_dp=1 1669 do ij=1,npaw_ij_in 1670 nspden=paw_ij_in(ij)%nspden 1671 buf_int(indx_int)=my_atmtab(ij) ;indx_int=indx_int+1 1672 buf_int(indx_int)=paw_ij_in(ij)%cplex_rf ;indx_int=indx_int+1 1673 buf_int(indx_int)=paw_ij_in(ij)%cplex_dij ;indx_int=indx_int+1 1674 buf_int(indx_int)=paw_ij_in(ij)%itypat ;indx_int=indx_int+1 1675 buf_int(indx_int)=paw_ij_in(ij)%nspden ;indx_int=indx_int+1 1676 buf_int(indx_int)=paw_ij_in(ij)%nsppol ;indx_int=indx_int+1 1677 buf_int(indx_int)=paw_ij_in(ij)%lmn_size ;indx_int=indx_int+1 1678 buf_int(indx_int)=paw_ij_in(ij)%lmn2_size ;indx_int=indx_int+1 1679 buf_int(indx_int)=paw_ij_in(ij)%ndij ;indx_int=indx_int+1 1680 buf_int(indx_int)=paw_ij_in(ij)%has_dij ;indx_int=indx_int+1 1681 buf_int(indx_int)=paw_ij_in(ij)%has_dij0 ;indx_int=indx_int+1 1682 buf_int(indx_int)=paw_ij_in(ij)%has_dijexxc ;indx_int=indx_int+1 1683 buf_int(indx_int)=paw_ij_in(ij)%has_dijfock ;indx_int=indx_int+1 1684 buf_int(indx_int)=paw_ij_in(ij)%has_dijfr ;indx_int=indx_int+1 1685 buf_int(indx_int)=paw_ij_in(ij)%has_dijhartree ;indx_int=indx_int+1 1686 buf_int(indx_int)=paw_ij_in(ij)%has_dijhat ;indx_int=indx_int+1 1687 buf_int(indx_int)=paw_ij_in(ij)%has_dijnd ;indx_int=indx_int+1 1688 buf_int(indx_int)=paw_ij_in(ij)%has_dijso ;indx_int=indx_int+1 1689 buf_int(indx_int)=paw_ij_in(ij)%has_dijU ;indx_int=indx_int+1 1690 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc ;indx_int=indx_int+1 1691 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc_hat ;indx_int=indx_int+1 1692 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc_val ;indx_int=indx_int+1 1693 buf_int(indx_int)=paw_ij_in(ij)%has_exexch_pot ;indx_int=indx_int+1 1694 buf_int(indx_int)=paw_ij_in(ij)%has_pawu_occ ;indx_int=indx_int+1 1695 lmn2_size=paw_ij_in(ij)%lmn2_size 1696 cplxdij_lmn2_size=paw_ij_in(ij)%cplex_dij*lmn2_size 1697 cplxrf_lmn2_size=paw_ij_in(ij)%cplex_rf*lmn2_size 1698 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij_in(ij)%cplex_rf 1699 ndij=paw_ij_in(ij)%ndij 1700 if (paw_ij_in(ij)%has_dij==2) then 1701 ii=cplxdijrf_lmn2_size*ndij 1702 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dij,(/ii/)) 1703 indx_dp=indx_dp+ii 1704 end if 1705 if (paw_ij_in(ij)%has_dij0==2) then 1706 ii=lmn2_size 1707 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij_in(ij)%dij0(1:ii) 1708 indx_dp=indx_dp+ii 1709 end if 1710 if (paw_ij_in(ij)%has_dijexxc==2) then 1711 ii=cplxdij_lmn2_size*ndij 1712 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijexxc,(/ii/)) 1713 indx_dp=indx_dp+ii 1714 end if 1715 if (paw_ij_in(ij)%has_dijfock==2) then 1716 ii=cplxdij_lmn2_size*ndij 1717 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijfock,(/ii/)) 1718 indx_dp=indx_dp+ii 1719 end if 1720 if (paw_ij_in(ij)%has_dijfr==2) then 1721 ii=cplxdijrf_lmn2_size*ndij 1722 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijfr,(/ii/)) 1723 indx_dp=indx_dp+ii 1724 end if 1725 if (paw_ij_in(ij)%has_dijhartree==2) then 1726 ii=cplxrf_lmn2_size 1727 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij_in(ij)%dijhartree(1:ii) 1728 indx_dp=indx_dp+ii 1729 end if 1730 if (paw_ij_in(ij)%has_dijhat==2) then 1731 ii=cplxdijrf_lmn2_size*ndij 1732 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijhat,(/ii/)) 1733 indx_dp=indx_dp+ii 1734 end if 1735 if (paw_ij_in(ij)%has_dijnd==2) then 1736 ii=cplxdij_lmn2_size*ndij 1737 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijnd,(/ii/)) 1738 indx_dp=indx_dp+ii 1739 end if 1740 if (paw_ij_in(ij)%has_dijso==2) then 1741 ii=cplxdijrf_lmn2_size*ndij 1742 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijso,(/ii/)) 1743 indx_dp=indx_dp+ii 1744 end if 1745 if (paw_ij_in(ij)%has_dijU==2) then 1746 ii=cplxdijrf_lmn2_size*ndij 1747 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijU,(/ii/)) 1748 indx_dp=indx_dp+ii 1749 end if 1750 if (paw_ij_in(ij)%has_dijxc==2) then 1751 ii=cplxdijrf_lmn2_size*ndij 1752 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc,(/ii/)) 1753 indx_dp=indx_dp+ii 1754 end if 1755 if (paw_ij_in(ij)%has_dijxc_hat==2) then 1756 ii=cplxdij_lmn2_size*ndij 1757 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc_hat,(/ii/)) 1758 indx_dp=indx_dp+ii 1759 end if 1760 if (paw_ij_in(ij)%has_dijxc_val==2) then 1761 ii=cplxdij_lmn2_size*ndij 1762 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc_val,(/ii/)) 1763 indx_dp=indx_dp+ii 1764 end if 1765 if (paw_ij_in(ij)%has_pawu_occ>=1)then 1766 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,1) ;indx_int=indx_int+1 1767 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,2) ;indx_int=indx_int+1 1768 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,3) ;indx_int=indx_int+1 1769 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,4) ;indx_int=indx_int+1 1770 nocc=paw_ij_in(ij)%ndij 1771 buf_dp(indx_dp:indx_dp+nocc-1)=paw_ij_in(ij)%nocctot(1:nocc) 1772 indx_dp=indx_dp+nocc 1773 nocc=size(paw_ij_in(ij)%noccmmp) 1774 buf_dp(indx_dp:indx_dp+nocc-1)=reshape(paw_ij_in(ij)%noccmmp,(/nocc/)) 1775 indx_dp=indx_dp+nocc 1776 end if 1777 if (paw_ij_in(ij)%has_exexch_pot>=1) then 1778 sz1=size(paw_ij_in(ij)%vpawx,1);sz2=size(paw_ij_in(ij)%vpawx,2) 1779 sz3=size(paw_ij_in(ij)%vpawx,3) 1780 buf_int(indx_int)=sz1; indx_int=indx_int+1 1781 buf_int(indx_int)=sz2; indx_int=indx_int+1 1782 buf_int(indx_int)=sz3; indx_int=indx_int+1 1783 ii=sz1*sz2*sz3 1784 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%vpawx,(/(ii)/)) 1785 indx_dp=indx_dp+ii 1786 end if 1787 end do 1788 1789 !Check 1790 indx_int=indx_int-1;indx_dp=indx_dp-1 1791 if ((indx_int/=buf_int_size).or.(indx_dp/=buf_dp_size)) then 1792 write(msg,*) 'Wrong buffer sizes: buf_int_size=',buf_int_size, & 1793 ' indx_int_size=',indx_int,' buf_dp_size=',buf_dp_size , ' indx_dp=', indx_dp 1794 MSG_BUG(msg) 1795 end if 1796 1797 !Communicate (1 gather for integers, 1 gather for reals) 1798 LIBPAW_ALLOCATE(count_int,(nproc_atom)) 1799 LIBPAW_ALLOCATE(displ_int,(nproc_atom)) 1800 LIBPAW_ALLOCATE(count_dp ,(nproc_atom)) 1801 LIBPAW_ALLOCATE(displ_dp ,(nproc_atom)) 1802 LIBPAW_ALLOCATE(count_tot,(2*nproc_atom)) 1803 bufsz(1)=buf_int_size; bufsz(2)=buf_dp_size 1804 call xmpi_allgather(bufsz,2,count_tot,comm_atom,ierr) 1805 do ij=1,nproc_atom 1806 count_int(ij)=count_tot(2*ij-1) 1807 count_dp (ij)=count_tot(2*ij) 1808 end do 1809 displ_int(1)=0;displ_dp(1)=0 1810 do ij=2,nproc_atom 1811 displ_int(ij)=displ_int(ij-1)+count_int(ij-1) 1812 displ_dp (ij)=displ_dp (ij-1)+count_dp (ij-1) 1813 end do 1814 buf_int_size_all=sum(count_int) 1815 buf_dp_size_all =sum(count_dp) 1816 LIBPAW_DEALLOCATE(count_tot) 1817 if (master==-1.or.me_atom==master) then 1818 LIBPAW_ALLOCATE(buf_int_all,(buf_int_size_all)) 1819 LIBPAW_ALLOCATE(buf_dp_all ,(buf_dp_size_all)) 1820 else 1821 LIBPAW_ALLOCATE(buf_int_all,(0)) 1822 LIBPAW_ALLOCATE(buf_dp_all ,(0)) 1823 end if 1824 if (master==-1) then 1825 call xmpi_allgatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,comm_atom,ierr) 1826 call xmpi_allgatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,comm_atom,ierr) 1827 else 1828 call xmpi_gatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,master,comm_atom,ierr) 1829 call xmpi_gatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,master,comm_atom,ierr) 1830 end if 1831 LIBPAW_DEALLOCATE(count_int) 1832 LIBPAW_DEALLOCATE(displ_int) 1833 LIBPAW_DEALLOCATE(count_dp) 1834 LIBPAW_DEALLOCATE(displ_dp) 1835 1836 !Retrieve data from output buffer 1837 if (master==-1.or.me_atom==master) then 1838 indx_int=1;indx_dp=1 1839 call paw_ij_free(paw_ij_gathered) 1840 call paw_ij_nullify(paw_ij_gathered) 1841 do ij=1,npaw_ij_out 1842 iat=buf_int_all(indx_int) ;indx_int=indx_int+1 1843 paw_ij_gathered(iat)%cplex_rf=buf_int_all(indx_int) ;indx_int=indx_int+1 1844 paw_ij_gathered(iat)%cplex_dij=buf_int_all(indx_int) ;indx_int=indx_int+1 1845 paw_ij_gathered(iat)%itypat=buf_int_all(indx_int) ;indx_int=indx_int+1 1846 paw_ij_gathered(iat)%nspden=buf_int_all(indx_int) ;indx_int=indx_int+1 1847 paw_ij_gathered(iat)%nsppol=buf_int_all(indx_int) ;indx_int=indx_int+1 1848 paw_ij_gathered(iat)%lmn_size=buf_int_all(indx_int) ;indx_int=indx_int+1 1849 paw_ij_gathered(iat)%lmn2_size=buf_int_all(indx_int) ;indx_int=indx_int+1 1850 paw_ij_gathered(iat)%ndij=buf_int_all(indx_int) ;indx_int=indx_int+1 1851 paw_ij_gathered(iat)%has_dij=buf_int_all(indx_int) ;indx_int=indx_int+1 1852 paw_ij_gathered(iat)%has_dij0=buf_int_all(indx_int) ;indx_int=indx_int+1 1853 paw_ij_gathered(iat)%has_dijexxc=buf_int_all(indx_int) ;indx_int=indx_int+1 1854 paw_ij_gathered(iat)%has_dijfock=buf_int_all(indx_int) ;indx_int=indx_int+1 1855 paw_ij_gathered(iat)%has_dijfr=buf_int_all(indx_int) ;indx_int=indx_int+1 1856 paw_ij_gathered(iat)%has_dijhartree=buf_int_all(indx_int) ;indx_int=indx_int+1 1857 paw_ij_gathered(iat)%has_dijhat=buf_int_all(indx_int) ;indx_int=indx_int+1 1858 paw_ij_gathered(iat)%has_dijnd=buf_int_all(indx_int) ;indx_int=indx_int+1 1859 paw_ij_gathered(iat)%has_dijso=buf_int_all(indx_int) ;indx_int=indx_int+1 1860 paw_ij_gathered(iat)%has_dijU=buf_int_all(indx_int) ;indx_int=indx_int+1 1861 paw_ij_gathered(iat)%has_dijxc=buf_int_all(indx_int) ;indx_int=indx_int+1 1862 paw_ij_gathered(iat)%has_dijxc_hat=buf_int_all(indx_int) ;indx_int=indx_int+1 1863 paw_ij_gathered(iat)%has_dijxc_val=buf_int_all(indx_int) ;indx_int=indx_int+1 1864 paw_ij_gathered(iat)%has_exexch_pot=buf_int_all(indx_int) ;indx_int=indx_int+1 1865 paw_ij_gathered(iat)%has_pawu_occ=buf_int_all(indx_int) ;indx_int=indx_int+1 1866 if (paw_ij_gathered(iat)%has_pawu_occ>=1) then 1867 nocc1=buf_int_all(indx_int) ;indx_int=indx_int+1 1868 nocc2=buf_int_all(indx_int) ;indx_int=indx_int+1 1869 nocc3=buf_int_all(indx_int) ;indx_int=indx_int+1 1870 nocc4=buf_int_all(indx_int) ;indx_int=indx_int+1 1871 else 1872 nocc1=0;nocc2=0;nocc3=0;nocc4=0 1873 end if 1874 lmn2_size=paw_ij_gathered(iat)%lmn2_size 1875 cplxdij_lmn2_size=paw_ij_gathered(iat)%cplex_dij*lmn2_size 1876 cplxrf_lmn2_size=paw_ij_gathered(iat)%cplex_rf*lmn2_size 1877 cplxdijrf_lmn2_size=cplxdij_lmn2_size*paw_ij_gathered(iat)%cplex_rf 1878 ndij=paw_ij_gathered(iat)%ndij 1879 1880 if (paw_ij_gathered(iat)%has_dij>=1) then 1881 ii=cplxdijrf_lmn2_size 1882 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij,(ii,ndij)) 1883 if (paw_ij_gathered(iat)%has_dij==2) then 1884 paw_ij_gathered(iat)%dij(:,:)= & 1885 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1886 indx_dp=indx_dp+ii*ndij 1887 end if 1888 end if 1889 if (paw_ij_gathered(iat)%has_dij0 >=1) then 1890 ii=lmn2_size 1891 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij0,(ii)) 1892 if (paw_ij_gathered(iat)%has_dij0==2) then 1893 paw_ij_gathered(iat)%dij0(:)=buf_dp_all(indx_dp:indx_dp+ii-1) 1894 indx_dp=indx_dp+ii 1895 end if 1896 end if 1897 if (paw_ij_gathered(iat)%has_dijexxc >=1) then 1898 ii=cplxdij_lmn2_size 1899 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijexxc,(ii,ndij)) 1900 if (paw_ij_gathered(iat)%has_dijexxc==2) then 1901 paw_ij_gathered(iat)%dijexxc(:,:)= & 1902 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1903 indx_dp=indx_dp+ii*ndij 1904 end if 1905 end if 1906 if (paw_ij_gathered(iat)%has_dijfock >=1) then 1907 ii=cplxdij_lmn2_size 1908 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfock,(ii,ndij)) 1909 if (paw_ij_gathered(iat)%has_dijfock==2) then 1910 paw_ij_gathered(iat)%dijfock(:,:)= & 1911 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1912 indx_dp=indx_dp+ii*ndij 1913 end if 1914 end if 1915 if (paw_ij_gathered(iat)%has_dijfr >=1) then 1916 ii=cplxdijrf_lmn2_size 1917 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfr,(ii,ndij)) 1918 if (paw_ij_gathered(iat)%has_dijfr==2) then 1919 paw_ij_gathered(iat)%dijfr(:,:)= & 1920 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1921 indx_dp=indx_dp+ii*ndij 1922 end if 1923 end if 1924 if (paw_ij_gathered(iat)%has_dijhartree >=1) then 1925 ii=cplxrf_lmn2_size 1926 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhartree,(ii)) 1927 if (paw_ij_gathered(iat)%has_dijhartree==2) then 1928 paw_ij_gathered(iat)%dijhartree(:)=buf_dp_all(indx_dp:indx_dp+ii-1) 1929 indx_dp=indx_dp+ii 1930 end if 1931 end if 1932 if (paw_ij_gathered(iat)%has_dijhat >=1) then 1933 ii=cplxdijrf_lmn2_size 1934 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhat,(ii,ndij)) 1935 if (paw_ij_gathered(iat)%has_dijhat==2) then 1936 paw_ij_gathered(iat)%dijhat(:,:)= & 1937 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1938 indx_dp=indx_dp+ii*ndij 1939 end if 1940 endif 1941 if (paw_ij_gathered(iat)%has_dijnd >=1) then 1942 ii=cplxdij_lmn2_size 1943 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijnd,(ii,ndij)) 1944 if (paw_ij_gathered(iat)%has_dijnd==2) then 1945 paw_ij_gathered(iat)%dijnd(:,:)= & 1946 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1947 indx_dp=indx_dp+ii*ndij 1948 end if 1949 end if 1950 if (paw_ij_gathered(iat)%has_dijso >=1) then 1951 ii=cplxdijrf_lmn2_size 1952 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijso,(ii,ndij)) 1953 if (paw_ij_gathered(iat)%has_dijso==2) then 1954 paw_ij_gathered(iat)%dijso(:,:)= & 1955 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1956 indx_dp=indx_dp+ii*ndij 1957 end if 1958 end if 1959 if (paw_ij_gathered(iat)%has_dijU >=1) then 1960 ii=cplxdijrf_lmn2_size 1961 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijU,(ii,ndij)) 1962 if (paw_ij_gathered(iat)%has_dijU==2) then 1963 paw_ij_gathered(iat)%dijU(:,:)= & 1964 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1965 indx_dp=indx_dp+ii*ndij 1966 end if 1967 end if 1968 if (paw_ij_gathered(iat)%has_dijxc >=1) then 1969 ii=cplxdijrf_lmn2_size 1970 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc,(ii,ndij)) 1971 if (paw_ij_gathered(iat)%has_dijxc==2) then 1972 paw_ij_gathered(iat)%dijxc(:,:)= & 1973 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1974 indx_dp=indx_dp+ii*ndij 1975 end if 1976 end if 1977 if (paw_ij_gathered(iat)%has_dijxc_hat >=1) then 1978 ii=cplxdij_lmn2_size 1979 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_hat,(ii,ndij)) 1980 if (paw_ij_gathered(iat)%has_dijxc_hat==2) then 1981 paw_ij_gathered(iat)%dijxc_hat(:,:)= & 1982 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1983 indx_dp=indx_dp+ii*ndij 1984 end if 1985 end if 1986 if (paw_ij_gathered(iat)%has_dijxc_val >=1) then 1987 ii=cplxdij_lmn2_size 1988 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_val,(ii,ndij)) 1989 if (paw_ij_gathered(iat)%has_dijxc_val==2) then 1990 paw_ij_gathered(iat)%dijxc_val(:,:)= & 1991 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1992 indx_dp=indx_dp+ii*ndij 1993 end if 1994 end if 1995 if (paw_ij_gathered(iat)%has_pawu_occ >=1) then 1996 nocc=paw_ij_gathered(iat)%ndij 1997 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%nocctot,(nocc)) 1998 paw_ij_gathered(iat)%nocctot(1:nocc)=buf_dp_all(indx_dp:indx_dp+nocc-1) 1999 indx_dp=indx_dp+nocc 2000 nocc=nocc1*nocc2*nocc3*nocc4 2001 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%noccmmp,(nocc1,nocc2,nocc3,nocc4)) 2002 paw_ij_gathered(iat)%noccmmp(1:nocc1,1:nocc2,1:nocc3,1:nocc4)= & 2003 & reshape(buf_dp_all(indx_dp:indx_dp+nocc-1),(/nocc1,nocc2,nocc3,nocc4/)) 2004 indx_dp=indx_dp+nocc 2005 end if 2006 if (paw_ij_gathered(iat)%has_exexch_pot >=1) then 2007 sz1=buf_int_all(indx_int);indx_int=indx_int+1 2008 sz2=buf_int_all(indx_int);indx_int=indx_int+1 2009 sz3=buf_int_all(indx_int);indx_int=indx_int+1 2010 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%vpawx,(sz1,sz2,sz3)) 2011 if (paw_ij_gathered(iat)%has_exexch_pot == 2) then 2012 paw_ij_gathered(iat)%vpawx(:,:,:)=& 2013 & reshape(buf_dp_all(indx_dp:indx_dp+sz1*sz2*sz3-1),(/sz1,sz2,sz3/)) 2014 indx_dp=indx_dp+sz1*sz2*sz3 2015 end if 2016 end if 2017 end do 2018 indx_int=indx_int-1;indx_dp=indx_dp-1 2019 if ((indx_int/=buf_int_size_all).or.(indx_dp/=buf_dp_size_all)) then 2020 write(msg,*) 'Wrong buffer sizes: buf_int_size_all=',buf_int_size_all, & 2021 & ' indx_int=',indx_int, ' buf_dp_size_all=',buf_dp_size_all,' indx_dp=',indx_dp 2022 MSG_BUG(msg) 2023 end if 2024 end if 2025 2026 !Free memory 2027 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 2028 LIBPAW_DEALLOCATE(buf_int) 2029 LIBPAW_DEALLOCATE(buf_dp) 2030 LIBPAW_DEALLOCATE(buf_int_all) 2031 LIBPAW_DEALLOCATE(buf_dp_all) 2032 2033 end subroutine paw_ij_gather
m_paw_ij/paw_ij_init [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_init
FUNCTION
Initialize a Paw_ij data type.
INPUTS
cplex=1 if no phase is applied (GS), 2 if a exp(-iqr) phase is applied (Response Function at q<>0) natom=Number of atoms. ntypat=Number of types of atoms in cell. nspinor=number of spinor components nsppol=Number of independent spin polarizations. nspden=Number of spin-density components pawspnorb=1 if spin-orbit coupling is activated typat(natom)=Type of each atom Pawtab(ntypat)<type(pawtab_type)>=PAW tabulated starting data OPTIONAL INPUTS has_dij=1 to allocate Paw_ij%dij, 0 otherwise (default) has_dij0=1 to allocate Paw_ij%dij0, 0 otherwise (default) has_dijfr=1 to allocate Paw_ij%dijfr, 0 otherwise (default) has_dijhat=1 to allocate Paw_ij%dijhat, 0 otherwise (default) has_dijxc=1 to allocate Paw_ij%dijxc, 0 otherwise (default) has_dijxc_hat=1 to allocate Paw_ij%dijxc_hat, 0 otherwise (default) has_dijxc_val=1 to allocate Paw_ij%dijxc_val, 0 otherwise (default) has_dijhartree=1 to allocate Paw_ij%dijhartree, 0 otherwise (default) has_dijfock=1 to allocate Paw_ij%dijfock, 0 otherwise (default) has_dijnd=1 to allocate Paw_ij%dijnd, used only if some nucdipmom /= 0; otherwise 0 (default) has_dijso=1 to allocate Paw_ij%dijso, used only if pawspnorb>0. 0 otherwise (default) has_dijU=1 to allocate Paw_ij%dijU, used only if Pawtab(itypat)%usepawu>0. 0 otherwise (default). has_dijexxc=to allocate Paw_ij%dijxx, 0 otherwise (default) has_exexch_pot=1 to allocate potential used in PAW+(local exact exchange) formalism, 0 otherwise (default) has_pawu_occ=1 to allocate occupations used in PAW+U formalism, 0 otherwise (default) nucdipmom(3,natom)= (optional) array of nuclear dipole moments at atomic sites mpi_atmtab(:)=indexes of the atoms treated by current proc comm_atom=MPI communicator over atoms
OUTPUT
Paw_ij(natom)<type(paw_ij_type)>=data structure containing PAW arrays given on (i,j) channels. In output all the basic dimensions are defined and the arrays are allocated according to the input variables.
PARENTS
bethe_salpeter,d2frnl,dfpt_nstpaw,dfpt_rhofermi,dfpt_scfcv,ldau_self m_energy,paw_qpscgw,respfn,scfcv,screening,sigma
CHILDREN
SOURCE
324 subroutine paw_ij_init(Paw_ij,cplex,nspinor,nsppol,nspden,pawspnorb,natom,ntypat,typat,Pawtab,& 325 & has_dij,has_dij0,has_dijfock,has_dijfr,has_dijhartree,has_dijhat,& ! Optional 326 & has_dijxc,has_dijxc_hat,has_dijxc_val,has_dijnd,has_dijso,has_dijU,has_dijexxc,& ! Optional 327 & has_exexch_pot,has_pawu_occ,nucdipmom,& ! Optional 328 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 329 330 331 !This section has been created automatically by the script Abilint (TD). 332 !Do not modify the following lines by hand. 333 #undef ABI_FUNC 334 #define ABI_FUNC 'paw_ij_init' 335 !End of the abilint section 336 337 implicit none 338 339 !Arguments ------------------------------------ 340 !scalars 341 integer,intent(in) :: cplex,nspinor,nspden,nsppol,natom,ntypat,pawspnorb 342 integer,optional,intent(in) :: has_dij,has_dij0,has_dijfr,has_dijhat,has_dijxc,has_dijxc_hat,has_dijxc_val 343 integer,optional,intent(in) :: has_dijnd,has_dijso,has_dijhartree,has_dijfock,has_dijU,has_dijexxc 344 integer,optional,intent(in) :: has_exexch_pot,has_pawu_occ 345 integer,optional,intent(in) :: comm_atom 346 347 !arrays 348 integer,intent(in) :: typat(natom) 349 integer,optional,target,intent(in) :: mpi_atmtab(:) 350 real(dp),optional,intent(in) :: nucdipmom(3,natom) 351 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 352 type(Pawtab_type),intent(in) :: Pawtab(ntypat) 353 354 !Local variables------------------------------- 355 !scalars 356 integer :: cplex_dij,cplex_rf,iat,iat_tot,itypat,lmn2_size,my_comm_atom,my_natom,ndij 357 logical :: my_atmtab_allocated,paral_atom,with_nucdipmom 358 !arrays 359 integer,pointer :: my_atmtab(:) 360 361 ! ************************************************************************* 362 363 !@Paw_ij_type 364 365 with_nucdipmom=.false.;if (present(nucdipmom)) with_nucdipmom=any(abs(nucdipmom)>tol8) 366 367 !Set up parallelism over atoms 368 my_natom=size(paw_ij);if (my_natom==0) return 369 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 370 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 371 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 372 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 373 374 do iat=1,my_natom 375 iat_tot=iat;if (paral_atom) iat_tot=my_atmtab(iat) 376 itypat=typat(iat_tot) 377 378 cplex_dij=1;if ((nspinor==2).or.(with_nucdipmom)) cplex_dij=2 379 cplex_rf=cplex 380 381 lmn2_size =Pawtab(itypat)%lmn2_size 382 Paw_ij(iat)%cplex_rf =cplex_rf 383 Paw_ij(iat)%cplex_dij =cplex_dij 384 Paw_ij(iat)%itypat =itypat 385 Paw_ij(iat)%nspden =nspden 386 Paw_ij(iat)%nsppol =nsppol 387 Paw_ij(iat)%lmn_size =Pawtab(itypat)%lmn_size 388 Paw_ij(iat)%lmn2_size =lmn2_size 389 Paw_ij(iat)%ndij =MAX(nspinor**2,nspden) 390 391 ndij=Paw_ij(iat)%ndij 392 393 ! ================================== 394 ! === Allocations (all optional) === 395 ! ================================== 396 397 ! === Allocation for total Dij === 398 Paw_ij(iat)%has_dij=0 399 if (PRESENT(has_dij)) then 400 if (has_dij/=0) then 401 Paw_ij(iat)%has_dij=1 402 LIBPAW_ALLOCATE(Paw_ij(iat)%dij,(cplex_rf*cplex_dij*lmn2_size,ndij)) 403 Paw_ij(iat)%dij(:,:)=zero 404 end if 405 end if 406 407 ! === Allocation for atomic Dij === 408 Paw_ij(iat)%has_dij0=0 409 if (PRESENT(has_dij0)) then 410 if (has_dij0/=0) then 411 Paw_ij(iat)%has_dij0=1 412 LIBPAW_ALLOCATE(Paw_ij(iat)%dij0,(lmn2_size)) 413 Paw_ij(iat)%dij0(:)=zero 414 end if 415 end if 416 417 ! === Allocation for Dij local exact exchange === 418 Paw_ij(iat)%has_dijexxc=0 419 if (PRESENT(has_dijexxc)) then 420 if (has_dijexxc/=0.and.Pawtab(itypat)%useexexch>0) then 421 Paw_ij(iat)%has_dijexxc=1 422 LIBPAW_ALLOCATE(Paw_ij(iat)%dijexxc,(cplex_dij*lmn2_size,ndij)) 423 Paw_ij(iat)%dijexxc(:,:)=zero 424 end if 425 end if 426 427 ! === Allocation for Dij_Fock === 428 Paw_ij(iat)%has_dijfock=0 429 if (PRESENT(has_dijfock)) then 430 if (has_dijfock/=0) then 431 Paw_ij(iat)%has_dijfock=1 432 LIBPAW_ALLOCATE(Paw_ij(iat)%dijfock,(cplex_dij*lmn2_size,ndij)) 433 Paw_ij(iat)%dijfock(:,:)=zero 434 end if 435 end if 436 437 ! === Allocation for frozen part of 1st-order Dij === 438 Paw_ij(iat)%has_dijfr=0 439 if (PRESENT(has_dijfr)) then 440 if (has_dijfr/=0) then 441 Paw_ij(iat)%has_dijfr=1 442 LIBPAW_ALLOCATE(Paw_ij(iat)%dijfr,(cplex_rf*cplex_dij*lmn2_size,ndij)) 443 Paw_ij(iat)%dijfr(:,:)=zero 444 end if 445 end if 446 447 ! === Allocation for Dij_Hartree === 448 Paw_ij(iat)%has_dijhartree=0 449 if (PRESENT(has_dijhartree)) then 450 if (has_dijhartree/=0) then 451 Paw_ij(iat)%has_dijhartree=1 452 LIBPAW_ALLOCATE(Paw_ij(iat)%dijhartree,(cplex_rf*lmn2_size)) 453 Paw_ij(iat)%dijhartree(:)=zero 454 end if 455 end if 456 457 ! === Allocation for Dij_hat === 458 Paw_ij(iat)%has_dijhat=0 459 if (PRESENT(has_dijhat)) then 460 if (has_dijhat/=0) then 461 Paw_ij(iat)%has_dijhat=1 462 LIBPAW_ALLOCATE(Paw_ij(iat)%dijhat,(cplex_rf*cplex_dij*lmn2_size,ndij)) 463 Paw_ij(iat)%dijhat(:,:)=zero 464 end if 465 end if 466 467 ! === Allocation for Dij nuclear dipole moment === 468 Paw_ij(iat)%has_dijnd=0 469 if (PRESENT(has_dijnd)) then 470 if (has_dijnd/=0.and.with_nucdipmom) then 471 Paw_ij(iat)%has_dijnd=1 472 LIBPAW_ALLOCATE(Paw_ij(iat)%dijnd,(cplex_dij*lmn2_size,ndij)) 473 Paw_ij(iat)%dijnd(:,:)=zero 474 end if 475 end if 476 477 ! === Allocation for Dij_SO === 478 Paw_ij(iat)%has_dijso=0 479 if (PRESENT(has_dijso)) then 480 if (has_dijso/=0.and.pawspnorb>0) then 481 Paw_ij(iat)%has_dijso=1 482 LIBPAW_ALLOCATE(Paw_ij(iat)%dijso,(cplex_rf*cplex_dij*lmn2_size,ndij)) 483 Paw_ij(iat)%dijso(:,:)=zero 484 end if 485 end if 486 487 ! === Allocation for Dij_U_val === 488 Paw_ij(iat)%has_dijU=0 489 if (PRESENT(has_dijU)) then 490 if (has_dijU/=0) then 491 Paw_ij(iat)%has_dijU=1 492 LIBPAW_ALLOCATE(Paw_ij(iat)%dijU,(cplex_rf*cplex_dij*lmn2_size,ndij)) 493 Paw_ij(iat)%dijU(:,:)=zero 494 end if 495 end if 496 497 ! === Allocation for total Dij_XC === 498 Paw_ij(iat)%has_dijxc=0 499 if (PRESENT(has_dijxc)) then 500 if (has_dijxc/=0) then 501 Paw_ij(iat)%has_dijxc=1 502 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc,(cplex_rf*cplex_dij*lmn2_size,ndij)) 503 Paw_ij(iat)%dijxc(:,:)=zero 504 end if 505 end if 506 507 ! === Allocation for total Dij_XC_hat === 508 Paw_ij(iat)%has_dijxc_hat=0 509 if (PRESENT(has_dijxc_hat)) then 510 if (has_dijxc_hat/=0) then 511 Paw_ij(iat)%has_dijxc_hat=1 512 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc_hat,(cplex_dij*lmn2_size,ndij)) 513 Paw_ij(iat)%dijxc_hat(:,:)=zero 514 end if 515 end if 516 517 ! === Allocation for total Dij_XC_val === 518 Paw_ij(iat)%has_dijxc_val=0 519 if (PRESENT(has_dijxc_val)) then 520 if (has_dijxc_val/=0) then 521 Paw_ij(iat)%has_dijxc_val=1 522 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc_val,(cplex_dij*lmn2_size,ndij)) 523 Paw_ij(iat)%dijxc_val(:,:)=zero 524 end if 525 end if 526 527 ! === Allocation for PAW+U occupations === 528 Paw_ij(iat)%has_pawu_occ=0 529 if (PRESENT(has_pawu_occ)) then 530 if (has_pawu_occ/=0.and.Pawtab(itypat)%usepawu>0) then 531 Paw_ij(iat)%has_pawu_occ=1 532 LIBPAW_ALLOCATE(Paw_ij(iat)%noccmmp,(cplex_dij,2*Pawtab(itypat)%lpawu+1,2*Pawtab(itypat)%lpawu+1,ndij)) 533 LIBPAW_ALLOCATE(Paw_ij(iat)%nocctot,(ndij)) 534 Paw_ij(iat)%noccmmp(:,:,:,:)=zero 535 Paw_ij(iat)%nocctot(:)=zero 536 end if 537 end if 538 539 ! === Allocation for PAW+LEXX potential === 540 Paw_ij(iat)%has_exexch_pot=0 541 if (PRESENT(has_exexch_pot)) then 542 if (has_exexch_pot/=0.and.Pawtab(itypat)%useexexch>0) then 543 Paw_ij(iat)%has_exexch_pot=1 544 ! TODO solve issue with first dimension 545 LIBPAW_ALLOCATE(Paw_ij(iat)%vpawx,(1,lmn2_size,nspden)) 546 Paw_ij(iat)%vpawx(:,:,:)=zero 547 end if 548 end if 549 550 end do 551 552 !Destroy atom table used for parallelism 553 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 554 555 end subroutine paw_ij_init
m_paw_ij/paw_ij_isendreceive_getbuffer [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_isendreceive_getbuffer
FUNCTION
Fill a paw_ij structure with the buffers received in a receive operation This buffer should have been first extracted by a call to paw_ij_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 npaw_ij_send= number of sent atoms
OUTPUT
paw_ij= output datastructure filled with buffers receive in a receive operation
PARENTS
m_paw_ij
CHILDREN
SOURCE
2522 subroutine paw_ij_isendreceive_getbuffer(paw_ij,npaw_ij_send,atm_indx_recv,buf_int,buf_dp) 2523 2524 2525 !This section has been created automatically by the script Abilint (TD). 2526 !Do not modify the following lines by hand. 2527 #undef ABI_FUNC 2528 #define ABI_FUNC 'paw_ij_isendreceive_getbuffer' 2529 !End of the abilint section 2530 2531 implicit none 2532 2533 !Arguments ------------------------------------ 2534 !scalars 2535 integer,intent(in) :: npaw_ij_send 2536 !arrays 2537 integer,intent(in):: atm_indx_recv(:),buf_int(:) 2538 real(dp),intent(in):: buf_dp(:) 2539 type(paw_ij_type),target,intent(inout) :: paw_ij(:) 2540 2541 !Local variables------------------------------- 2542 !scalars 2543 integer :: buf_dp_size,buf_int_size 2544 integer :: cplxrf_lmn2_size,cplxdij_lmn2_size,cplxdijrf_lmn2_size 2545 integer :: iat,iatom_tot,ii,ij,indx_dp,indx_int,ndij,nocc,nocc1,nocc2,nocc3,nocc4 2546 integer :: lmn2_size,sz1,sz2,sz3 2547 character(len=500) :: msg 2548 type(Paw_ij_type),pointer :: paw_ij1 2549 !arrays 2550 2551 ! ********************************************************************* 2552 2553 buf_int_size=size(buf_int) 2554 buf_dp_size=size(buf_dp) 2555 indx_int=1;indx_dp=1 2556 2557 do ij=1,npaw_ij_send 2558 iatom_tot=buf_int(indx_int) ;indx_int=indx_int+1 2559 iat= atm_indx_recv(iatom_tot) 2560 paw_ij1=>paw_ij(iat) 2561 paw_ij1%cplex_rf=buf_int(indx_int) ;indx_int=indx_int+1 2562 paw_ij1%cplex_dij=buf_int(indx_int) ;indx_int=indx_int+1 2563 paw_ij1%itypat=buf_int(indx_int) ;indx_int=indx_int+1 2564 paw_ij1%nspden=buf_int(indx_int) ;indx_int=indx_int+1 2565 paw_ij1%nsppol=buf_int(indx_int) ;indx_int=indx_int+1 2566 paw_ij1%lmn_size=buf_int(indx_int) ;indx_int=indx_int+1 2567 paw_ij1%lmn2_size=buf_int(indx_int) ;indx_int=indx_int+1 2568 paw_ij1%ndij=buf_int(indx_int) ;indx_int=indx_int+1 2569 paw_ij1%has_dij=buf_int(indx_int) ;indx_int=indx_int+1 2570 paw_ij1%has_dij0=buf_int(indx_int) ;indx_int=indx_int+1 2571 paw_ij1%has_dijexxc=buf_int(indx_int) ;indx_int=indx_int+1 2572 paw_ij1%has_dijfock=buf_int(indx_int) ;indx_int=indx_int+1 2573 paw_ij1%has_dijfr=buf_int(indx_int) ;indx_int=indx_int+1 2574 paw_ij1%has_dijhartree=buf_int(indx_int) ;indx_int=indx_int+1 2575 paw_ij1%has_dijhat=buf_int(indx_int) ;indx_int=indx_int+1 2576 paw_ij1%has_dijnd=buf_int(indx_int) ;indx_int=indx_int+1 2577 paw_ij1%has_dijso=buf_int(indx_int) ;indx_int=indx_int+1 2578 paw_ij1%has_dijU=buf_int(indx_int) ;indx_int=indx_int+1 2579 paw_ij1%has_dijxc=buf_int(indx_int) ;indx_int=indx_int+1 2580 paw_ij1%has_dijxc_hat=buf_int(indx_int) ;indx_int=indx_int+1 2581 paw_ij1%has_dijxc_val=buf_int(indx_int) ;indx_int=indx_int+1 2582 paw_ij1%has_exexch_pot=buf_int(indx_int) ;indx_int=indx_int+1 2583 paw_ij1%has_pawu_occ=buf_int(indx_int) ;indx_int=indx_int+1 2584 if (paw_ij1%has_pawu_occ>=1) then 2585 nocc1=buf_int(indx_int) ;indx_int=indx_int+1 2586 nocc2=buf_int(indx_int) ;indx_int=indx_int+1 2587 nocc3=buf_int(indx_int) ;indx_int=indx_int+1 2588 nocc4=buf_int(indx_int) ;indx_int=indx_int+1 2589 else 2590 nocc1=0;nocc2=0;nocc3=0;nocc4=0 2591 end if 2592 lmn2_size=paw_ij1%lmn2_size 2593 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2594 cplxrf_lmn2_size=paw_ij1%cplex_rf*lmn2_size 2595 cplxdijrf_lmn2_size=cplxrf_lmn2_size*paw_ij1%cplex_rf 2596 ndij=paw_ij1%ndij 2597 2598 if (paw_ij1%has_dij>=1) then 2599 ii=cplxdijrf_lmn2_size 2600 LIBPAW_ALLOCATE(paw_ij1%dij,(ii,ndij)) 2601 if (paw_ij1%has_dij==2) then 2602 paw_ij1%dij(:,:)= & 2603 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2604 indx_dp=indx_dp+ii*ndij 2605 end if 2606 end if 2607 if (paw_ij1%has_dij0 >=1) then 2608 ii=lmn2_size 2609 LIBPAW_ALLOCATE(paw_ij1%dij0,(ii)) 2610 if (paw_ij1%has_dij0==2) then 2611 paw_ij1%dij0(:)=buf_dp(indx_dp:indx_dp+ii-1) 2612 indx_dp=indx_dp+ii 2613 end if 2614 end if 2615 if (paw_ij1%has_dijexxc >=1) then 2616 ii=cplxdij_lmn2_size 2617 LIBPAW_ALLOCATE(paw_ij1%dijexxc,(ii,ndij)) 2618 if (paw_ij1%has_dijexxc==2) then 2619 paw_ij1%dijexxc(:,:)= & 2620 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2621 indx_dp=indx_dp+ii*ndij 2622 end if 2623 end if 2624 if (paw_ij1%has_dijfock >=1) then 2625 ii=cplxdij_lmn2_size 2626 LIBPAW_ALLOCATE(paw_ij1%dijfock,(ii,ndij)) 2627 if (paw_ij1%has_dijfock==2) then 2628 paw_ij1%dijfock(:,:)= & 2629 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2630 indx_dp=indx_dp+ii*ndij 2631 end if 2632 end if 2633 if (paw_ij1%has_dijfr >=1) then 2634 ii=cplxdijrf_lmn2_size 2635 LIBPAW_ALLOCATE(paw_ij1%dijfr,(ii,ndij)) 2636 if (paw_ij1%has_dijfr==2) then 2637 paw_ij1%dijfr(:,:)= & 2638 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2639 indx_dp=indx_dp+ii*ndij 2640 end if 2641 end if 2642 if (paw_ij1%has_dijhartree >=1) then 2643 ii=cplxrf_lmn2_size 2644 LIBPAW_ALLOCATE(paw_ij1%dijhartree,(ii)) 2645 if (paw_ij1%has_dijhartree==2) then 2646 paw_ij1%dijhartree(:)=buf_dp(indx_dp:indx_dp+ii-1) 2647 indx_dp=indx_dp+ii 2648 end if 2649 end if 2650 if (paw_ij1%has_dijhat >=1) then 2651 ii=cplxdijrf_lmn2_size 2652 LIBPAW_ALLOCATE(paw_ij1%dijhat,(ii,ndij)) 2653 if (paw_ij1%has_dijhat==2) then 2654 paw_ij1%dijhat(:,:)= & 2655 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2656 indx_dp=indx_dp+ii*ndij 2657 end if 2658 end if 2659 if (paw_ij1%has_dijnd >=1) then 2660 ii=cplxdij_lmn2_size 2661 LIBPAW_ALLOCATE(paw_ij1%dijnd,(ii,ndij)) 2662 if (paw_ij1%has_dijnd==2) then 2663 paw_ij1%dijnd(:,:)= & 2664 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2665 indx_dp=indx_dp+ii*ndij 2666 end if 2667 end if 2668 if (paw_ij1%has_dijso >=1) then 2669 ii=cplxdijrf_lmn2_size 2670 LIBPAW_ALLOCATE(paw_ij1%dijso,(ii,ndij)) 2671 if (paw_ij1%has_dijso==2) then 2672 paw_ij1%dijso(:,:)= & 2673 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2674 indx_dp=indx_dp+ii*ndij 2675 end if 2676 end if 2677 if (paw_ij1%has_dijU >=1) then 2678 ii=cplxdijrf_lmn2_size 2679 LIBPAW_ALLOCATE(paw_ij1%dijU,(ii,ndij)) 2680 if (paw_ij1%has_dijU==2) then 2681 paw_ij1%dijU(:,:)= & 2682 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2683 indx_dp=indx_dp+ii*ndij 2684 end if 2685 end if 2686 if (paw_ij1%has_dijxc >=1) then 2687 ii=cplxdijrf_lmn2_size 2688 LIBPAW_ALLOCATE(paw_ij1%dijxc,(ii,ndij)) 2689 if (paw_ij1%has_dijxc==2) then 2690 paw_ij1%dijxc(:,:)= & 2691 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2692 indx_dp=indx_dp+ii*ndij 2693 end if 2694 end if 2695 if (paw_ij1%has_dijxc_hat >=1) then 2696 ii=cplxdij_lmn2_size 2697 LIBPAW_ALLOCATE(paw_ij1%dijxc_hat,(ii,ndij)) 2698 if (paw_ij1%has_dijxc_hat==2) then 2699 paw_ij1%dijxc_hat(:,:)= & 2700 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2701 indx_dp=indx_dp+ii*ndij 2702 end if 2703 end if 2704 if (paw_ij1%has_dijxc_val >=1) then 2705 ii=cplxdij_lmn2_size 2706 LIBPAW_ALLOCATE(paw_ij1%dijxc_val,(ii,ndij)) 2707 if (paw_ij1%has_dijxc_val==2) then 2708 paw_ij1%dijxc_val(:,:)= & 2709 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2710 indx_dp=indx_dp+ii*ndij 2711 end if 2712 end if 2713 if (paw_ij1%has_pawu_occ >=1) then 2714 nocc=paw_ij1%ndij 2715 LIBPAW_ALLOCATE(paw_ij1%nocctot,(nocc)) 2716 paw_ij1%nocctot(1:nocc)=buf_dp(indx_dp:indx_dp+nocc-1) 2717 indx_dp=indx_dp+nocc 2718 nocc=nocc1*nocc2*nocc3*nocc4 2719 LIBPAW_ALLOCATE(paw_ij1%noccmmp,(nocc1,nocc2,nocc3,nocc4)) 2720 paw_ij1%noccmmp(1:nocc1,1:nocc2,1:nocc3,1:nocc4)= & 2721 & reshape(buf_dp(indx_dp:indx_dp+nocc-1),(/nocc1,nocc2,nocc3,nocc4/)) 2722 indx_dp=indx_dp+nocc 2723 end if 2724 if (paw_ij1%has_exexch_pot >=1) then 2725 sz1=buf_int(indx_int);indx_int=indx_int+1 2726 sz2=buf_int(indx_int);indx_int=indx_int+1 2727 sz3=buf_int(indx_int);indx_int=indx_int+1 2728 LIBPAW_ALLOCATE(paw_ij1%vpawx,(sz1,sz2,sz3)) 2729 if (paw_ij1%has_exexch_pot == 2) then 2730 paw_ij1%vpawx(:,:,:)=& 2731 & reshape(buf_dp(indx_dp:indx_dp+sz1*sz2*sz3-1),(/sz1,sz2,sz3/)) 2732 indx_dp=indx_dp+sz1*sz2*sz3 2733 end if 2734 end if 2735 end do 2736 2737 if ((indx_int/=1+buf_int_size).or.(indx_dp /=1+buf_dp_size)) then 2738 write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size 2739 MSG_BUG(msg) 2740 end if 2741 2742 end subroutine paw_ij_isendreceive_getbuffer
m_paw_ij/paw_ij_nullify [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_nullify
FUNCTION
Reset all flags in a paw_ij structure
SIDE EFFECTS
Paw_ij(:)<type(paw_ij_type)>=PAW arrays given on (i,j) channels.
PARENTS
bethe_salpeter,d2frnl,dfpt_nstpaw,dfpt_rhofermi,dfpt_scfcv,ldau_self m_energy,m_paw_ij,paw_qpscgw,pawprt,respfn,scfcv,screening,sigma
CHILDREN
SOURCE
695 subroutine paw_ij_nullify(Paw_ij) 696 697 698 !This section has been created automatically by the script Abilint (TD). 699 !Do not modify the following lines by hand. 700 #undef ABI_FUNC 701 #define ABI_FUNC 'paw_ij_nullify' 702 !End of the abilint section 703 704 implicit none 705 706 !Arguments ------------------------------------ 707 !arrays 708 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 709 710 !Local variables------------------------------- 711 integer :: iat,natom 712 713 ! ************************************************************************* 714 715 !@Paw_ij_type 716 717 ! MGPAW: This one could be removed/renamed, 718 ! variables can be initialized in the datatype declaration 719 ! Do we need to expose this in the public API? 720 721 natom=SIZE(Paw_ij(:));if (natom==0) return 722 723 ! Set all has_* flags to zero. 724 do iat=1,natom 725 ! === Set all has_* flags to zero === 726 Paw_ij(iat)%has_dij =0 727 Paw_ij(iat)%has_dij0 =0 728 Paw_ij(iat)%has_dijexxc =0 729 Paw_ij(iat)%has_dijfock =0 730 Paw_ij(iat)%has_dijfr =0 731 Paw_ij(iat)%has_dijhartree=0 732 Paw_ij(iat)%has_dijhat =0 733 Paw_ij(iat)%has_dijnd =0 734 Paw_ij(iat)%has_dijso =0 735 Paw_ij(iat)%has_dijU =0 736 Paw_ij(iat)%has_dijxc =0 737 Paw_ij(iat)%has_dijxc_hat =0 738 Paw_ij(iat)%has_dijxc_val =0 739 Paw_ij(iat)%has_exexch_pot=0 740 Paw_ij(iat)%has_pawu_occ =0 741 end do !iat 742 743 end subroutine paw_ij_nullify
m_paw_ij/paw_ij_print [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_print
FUNCTION
Print out the content of a paw_ij datastructure (Dij only)
INPUTS
[enunit]=governs the units to be used for the output of total Dij (0:Ha, 1:Ha+eV) [ipert]=only for DFPT: index of the perturbation (0 for ground state) [unit]=the unit number for output [pawprtvol]=verbosity level [pawspnorb]=1 if spin-orbit coupling is activated [mode_paral]=either "COLL" or "PERS" [mpi_atmtab(:)]=indexes of the atoms treated by current proc (can be computed here) [comm_atom]=MPI communicator over atoms (needed if parallelism over atoms is activated) [natom]=total number of atom (needed if parallelism over atoms is activated) if Paw_ij is distributed, natom is different from size(Paw_ij).
OUTPUT
(Only writing)
NOTES
PARENTS
m_pawdij,sigma
CHILDREN
SOURCE
1018 subroutine paw_ij_print(Paw_ij,unit,pawprtvol,pawspnorb,mode_paral,enunit,ipert, & 1019 & mpi_atmtab,comm_atom,natom) 1020 1021 1022 !This section has been created automatically by the script Abilint (TD). 1023 !Do not modify the following lines by hand. 1024 #undef ABI_FUNC 1025 #define ABI_FUNC 'paw_ij_print' 1026 !End of the abilint section 1027 1028 implicit none 1029 1030 !Arguments ------------------------------------ 1031 !scalars 1032 integer,optional,intent(in) :: enunit,ipert 1033 integer,optional,intent(in) :: comm_atom,natom 1034 integer,optional,intent(in) :: pawprtvol,pawspnorb 1035 integer,optional,intent(in) :: unit 1036 character(len=4),optional,intent(in) :: mode_paral 1037 !arrays 1038 integer,optional,target,intent(in) :: mpi_atmtab(:) 1039 type(Paw_ij_type),target,intent(in) :: Paw_ij(:) 1040 1041 !Local variables------------------------------- 1042 character(len=7),parameter :: dspin(6)=(/"up ","down ","up-up ","dwn-dwn","up-dwn ","dwn-up "/) 1043 !scalars 1044 integer :: cplex_dij,cplex_rf,iatom,iatom_tot,idij,idij_sym,klmn,lmn2_size,lmn_size,my_comm_atom,my_natom 1045 integer :: nspden,nsploop,nsppol,my_unt,ndij,tmp_cplex_dij,my_ipert,my_enunit,my_prtvol,size_paw_ij 1046 logical :: my_atmtab_allocated,paral_atom 1047 character(len=4) :: my_mode 1048 character(len=2000) :: msg 1049 !arrays 1050 integer :: idum(0) 1051 integer,pointer :: my_atmtab(:) 1052 real(dp),allocatable,target :: dij(:),dijs(:),dijh(:,:) 1053 real(dp),pointer :: dij2p(:),dij2p_(:) 1054 1055 ! ************************************************************************* 1056 1057 !@Paw_ij_type 1058 1059 size_paw_ij=SIZE(Paw_ij);if (size_paw_ij==0) return 1060 1061 my_unt =std_out ; if (PRESENT(unit )) my_unt =unit 1062 my_prtvol=0 ; if (PRESENT(pawprtvol )) my_prtvol=pawprtvol 1063 my_mode ='COLL' ; if (PRESENT(mode_paral)) my_mode =mode_paral 1064 my_ipert =0 ; if (PRESENT(ipert)) my_ipert =ipert 1065 my_enunit=0 ; if (PRESENT(enunit)) my_enunit=enunit 1066 my_natom=size_paw_ij; if (PRESENT(natom)) my_natom=natom 1067 1068 !Set up parallelism over atoms 1069 paral_atom=(present(comm_atom).and.my_natom/=size_paw_ij) 1070 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1071 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 1072 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,my_natom,my_natom_ref=size_paw_ij) 1073 1074 if (abs(my_prtvol)>=1) then 1075 if (my_ipert==0) then 1076 write(msg,'(4a)')ch10,' ==== Values of psp strength Dij (Hartree) ============' 1077 else 1078 write(msg,'(4a)')ch10,' ==== Values of psp strength Dij(1) (Hartree) =========' 1079 end if 1080 call wrtout(my_unt,msg,my_mode) 1081 end if 1082 1083 nsppol = Paw_ij(1)%nsppol 1084 nspden = Paw_ij(1)%nspden 1085 nsploop= nsppol; if (Paw_ij(1)%ndij==4) nsploop=4 1086 1087 do iatom=1,size_paw_ij 1088 1089 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1090 1091 lmn_size = Paw_ij(iatom)%lmn_size 1092 lmn2_size = Paw_ij(iatom)%lmn2_size 1093 cplex_dij = Paw_ij(iatom)%cplex_dij 1094 cplex_rf = Paw_ij(iatom)%cplex_rf 1095 ndij = Paw_ij(iatom)%ndij 1096 1097 ! ==================================== 1098 ! === Loop over density components === 1099 ! ==================================== 1100 do idij=1,nsploop 1101 1102 idij_sym=idij;if (ndij==4.and.idij>2) idij_sym=7-idij 1103 if (cplex_rf==2) then 1104 LIBPAW_ALLOCATE(dij,(2*lmn2_size)) 1105 LIBPAW_ALLOCATE(dijs,(2*lmn2_size)) 1106 end if 1107 1108 ! =================== Detailed output ===================================== 1109 if (ABS(my_prtvol)>=1.and.(iatom_tot==1.or.iatom_tot==my_natom.or.my_prtvol<0)) then 1110 1111 !Title 1112 if (nspden==2.and.nsppol==1) then 1113 write(msg,'(2a,i3,3a)')ch10,& 1114 & ' >>>>>>>>>> Atom ',iatom_tot,':',ch10,& 1115 & ' (antiferromagnetism case: only one spin component)' 1116 else if (paw_ij(iatom)%ndij==1) then 1117 write(msg, '(2a,i3,a)') ch10,& 1118 & ' >>>>>>>>>> Atom ',iatom_tot,':' 1119 else 1120 write(msg,'(2a,i3,3a)') ch10,& 1121 & ' >>>>>>>>>> Atom ',iatom_tot,' (component ',TRIM(dspin(idij+2*(nsploop/4))),'):' 1122 end if 1123 call wrtout(my_unt,msg,my_mode) 1124 1125 !Dij atomic 1126 if (Paw_ij(iatom)%has_dij0/=0.and.idij<=2.and.my_ipert<=0) then 1127 write(msg,'(a)') ' ************ Dij atomic (Dij0) ***********' 1128 call wrtout(my_unt,msg,my_mode) 1129 call pawio_print_ij(my_unt,Paw_ij(iatom)%dij0,lmn2_size,1,lmn_size,-1,idum,0,my_prtvol,idum,-1.d0,1,& 1130 & opt_sym=2,mode_paral=my_mode) 1131 end if 1132 1133 !Dij Local Exact Exchange 1134 if (Paw_ij(iatom)%has_dijexxc/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1135 write(msg,'(a)') ' ************* Dij_Local Exact exchange **********' 1136 call wrtout(my_unt,msg,my_mode) 1137 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijexxc) 1138 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1139 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1140 end if 1141 1142 !Dij Fock 1143 if (Paw_ij(iatom)%has_dijfock/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1144 write(msg,'(a)') ' ************* Dij_Fock **********' 1145 call wrtout(my_unt,msg,my_mode) 1146 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijfock) 1147 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1148 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1149 end if 1150 1151 !Dij Frozen (RF) 1152 if (Paw_ij(iatom)%has_dijfr/=0.and.(idij<=2.or.nspden==4).and.my_ipert>0) then 1153 write(msg,'(a)') ' ************** Dij(1) Frozen **************' 1154 call wrtout(my_unt,msg,my_mode) 1155 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dijfr) 1156 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1157 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1158 end if 1159 1160 !Dij Hartree 1161 if (Paw_ij(iatom)%has_dijhartree/=0.and.idij<=2) then 1162 if (my_ipert==0) then 1163 write(msg,'(a)') ' ************** Dij Hartree ***************' 1164 else 1165 write(msg,'(a)') ' ************* Dij(1) Hartree *************' 1166 end if 1167 call wrtout(my_unt,msg,my_mode) 1168 LIBPAW_ALLOCATE(dijh,(cplex_rf*lmn2_size,1)) 1169 dijh(:,1)=Paw_ij(iatom)%dijhartree(:) 1170 call get_dij_parts(1,cplex_rf,dijh) 1171 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0, & 1172 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1173 LIBPAW_DEALLOCATE(dijh) 1174 end if 1175 1176 !Dij Hat 1177 if (Paw_ij(iatom)%has_dijhat/=0.and.(idij<=2.or.nspden==4)) then 1178 if (my_ipert==0) then 1179 write(msg,'(a)') ' **************** Dij_hat *****************' 1180 else 1181 write(msg,'(a)') ' ***** Dij_hat(1) (incl. frozen Dij) ******' 1182 end if 1183 call wrtout(my_unt,msg,my_mode) 1184 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dijhat) 1185 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1186 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1187 end if 1188 1189 !Dij nuclear dipole 1190 if (Paw_ij(iatom)%has_dijnd/=0.and.my_ipert<=0) then 1191 write(msg,'(a)') ' *********** Dij Nuclear Dipole **********' 1192 call wrtout(my_unt,msg,my_mode) 1193 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijnd,always_img=.true.) 1194 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1195 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1196 end if 1197 1198 !Dij spin-orbit 1199 if (Paw_ij(iatom)%has_dijso/=0.and.my_ipert<=0) then 1200 write(msg,'(a)') ' ************** Dij SpinOrbit ************' 1201 call wrtout(my_unt,msg,my_mode) 1202 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dijso,always_img=.true.) 1203 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1204 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1205 end if 1206 1207 !Dij LDA+U 1208 if (Paw_ij(iatom)%has_dijU/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1209 write(msg,'(a)') ' ************* Dij_LDA+U (dijpawu) **********' 1210 call wrtout(my_unt,msg,my_mode) 1211 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%diju) 1212 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1213 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1214 end if 1215 1216 !Dij XC 1217 if (Paw_ij(iatom)%has_dijxc/=0.and.(idij<=2.or.nspden==4)) then 1218 if (my_ipert<=0) then 1219 write(msg,'(a)') ' ***************** Dij_xc *****************' 1220 else 1221 write(msg,'(a)') ' **************** Dij(1)_xc ***************' 1222 end if 1223 call wrtout(my_unt,msg,my_mode) 1224 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dijxc) 1225 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1226 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1227 end if 1228 1229 !Dij hat XC 1230 if (Paw_ij(iatom)%has_dijxc_hat/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1231 if (my_ipert<=0) then 1232 write(msg,'(a)') ' *************** Dijhat_xc ****************' 1233 else 1234 write(msg,'(a)') ' ************** Dij(1)hat_xc **************' 1235 end if 1236 call wrtout(my_unt,msg,my_mode) 1237 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijxc_hat) 1238 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1239 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1240 end if 1241 1242 !Dij XC val 1243 if (Paw_ij(iatom)%has_dijxc_val/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1244 write(msg, '(a)') ' *************** Dij_xc_val ***************' 1245 call wrtout(my_unt,msg,my_mode) 1246 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijxc_val) 1247 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1248 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1249 end if 1250 1251 !Dij TOTAL 1252 if (Paw_ij(iatom)%has_dij/=0) then 1253 if (my_ipert<=0) then 1254 write(msg,'(a)') ' ********** TOTAL Dij in Ha **********' 1255 else 1256 write(msg,'(a)') ' ********** TOTAL Dij(1) in Ha **********' 1257 end if 1258 call wrtout(my_unt,msg,my_mode) 1259 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dij,always_img=.true.) 1260 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1261 & my_prtvol,idum,50.d0*dble(3-2*idij),1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1262 if (my_enunit>0) then 1263 if (my_ipert<=0) then 1264 write(msg,'(a)') ' ********** TOTAL Dij in eV **********' 1265 else 1266 write(msg,'(a)') ' ********** TOTAL Dij(1) in eV **********' 1267 end if 1268 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1269 & my_prtvol,idum,-1._dp,2,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1270 end if 1271 end if 1272 1273 if (allocated(dij)) then 1274 LIBPAW_DEALLOCATE(dij) 1275 end if 1276 if (allocated(dijs)) then 1277 LIBPAW_DEALLOCATE(dijs) 1278 end if 1279 1280 end if !(ABS(my_prtvol)>=1.and.(iatom_tot==1.or.iatom_tot==my_natom.or.my_prtvol<0) 1281 1282 ! =================== Standard output ===================================== 1283 if ((abs(my_prtvol)==0).and.(iatom_tot==1.or.iatom_tot==my_natom)) then 1284 1285 !Title 1286 if (idij==1) then 1287 if (my_ipert<=0) then 1288 write(msg, '(2a,i6,a)') ch10,' ****** Psp strength Dij in Ha (atom ',iatom_tot,') *****' 1289 else 1290 write(msg, '(2a,i6,a)') ch10,' **** Psp strength Dij(1) in Ha (atom ',iatom_tot,') *****' 1291 end if 1292 if (nspden==2.and.nsppol==1) then 1293 write(msg,'(4a)') trim(msg),') *****',ch10,' (antiferromagnetism case: only one spin component)' 1294 end if 1295 call wrtout(my_unt,msg,my_mode) 1296 end if 1297 if (paw_ij(iatom)%ndij/=1) then 1298 write(msg,'(3a)') ' Component ',trim(dspin(idij+2*(nsploop/4))),':' 1299 call wrtout(my_unt,msg,my_mode) 1300 end if 1301 1302 !Dij TOTAL 1303 if (Paw_ij(iatom)%has_dij/=0) then 1304 call get_dij_parts(cplex_dij,cplex_rf,Paw_ij(iatom)%dij,always_img=.true.) 1305 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1306 & my_prtvol,idum,50.d0*dble(3-2*idij),1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1307 end if 1308 1309 end if 1310 1311 ! =================== End main loops ===================================== 1312 1313 if (allocated(dij)) then 1314 LIBPAW_DEALLOCATE(dij) 1315 end if 1316 if (allocated(dijs)) then 1317 LIBPAW_DEALLOCATE(dijs) 1318 end if 1319 1320 end do !idij 1321 end do !iat 1322 1323 call wrtout(my_unt,' ',my_mode) 1324 1325 !Small helper function 1326 contains 1327 1328 !Real and imaginary parts of phase. 1329 subroutine get_dij_parts(my_cplex_dij,my_cplex_rf,my_dij,always_img) 1330 1331 1332 !This section has been created automatically by the script Abilint (TD). 1333 !Do not modify the following lines by hand. 1334 #undef ABI_FUNC 1335 #define ABI_FUNC 'get_dij_parts' 1336 !End of the abilint section 1337 1338 integer,intent(in) :: my_cplex_dij,my_cplex_rf 1339 logical,intent(in),optional :: always_img 1340 real(dp),intent(in),target :: my_dij(:,:) 1341 integer :: my_idij,my_idij_sym,kk 1342 logical :: always_img_ 1343 always_img_=.false.;if(present(always_img)) always_img_=always_img 1344 my_idij=min(size(my_dij,2),idij) 1345 my_idij_sym=min(size(my_dij,2),idij_sym) 1346 if (my_cplex_rf==1) then 1347 if ((idij<=nsppol.or.idij==2).and.(.not.always_img_))then 1348 tmp_cplex_dij=1 1349 dij2p => my_dij(1:my_cplex_dij*lmn2_size:my_cplex_dij,my_idij) 1350 dij2p_ => dij2p 1351 else 1352 tmp_cplex_dij=my_cplex_dij 1353 dij2p => my_dij(1:my_cplex_dij*lmn2_size:1,my_idij) 1354 dij2p_ => my_dij(1:my_cplex_dij*lmn2_size:1,my_idij_sym) 1355 end if 1356 else 1357 tmp_cplex_dij=2 1358 if (my_cplex_dij==1) then 1359 do kk=1,lmn2_size 1360 dij(2*kk-1)= my_dij(kk,my_idij) 1361 dij(2*kk )= my_dij(kk+lmn2_size,my_idij) 1362 dijs(2*kk-1)= my_dij(kk,my_idij) 1363 dijs(2*kk )=-my_dij(kk+lmn2_size,my_idij) 1364 end do 1365 else 1366 do kk=1,lmn2_size 1367 dij(2*kk-1)= my_dij(2*kk-1,idij)-my_dij(2*kk +2*lmn2_size,my_idij) 1368 dij(2*kk )= my_dij(2*kk ,idij)+my_dij(2*kk-1+2*lmn2_size,my_idij) 1369 dijs(2*kk-1)= my_dij(2*kk-1,idij_sym)+my_dij(2*kk +2*lmn2_size,my_idij_sym) 1370 dijs(2*kk )= my_dij(2*kk ,idij_sym)-my_dij(2*kk-1+2*lmn2_size,my_idij_sym) 1371 end do 1372 end if 1373 dij2p => dij ; dij2p_ => dijs 1374 end if 1375 end subroutine get_dij_parts 1376 1377 end subroutine paw_ij_print
m_paw_ij/paw_ij_redistribute [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_redistribute
FUNCTION
Redistribute an array of paw_ij datastructures Input paw_ij is given on a MPI communicator Output paw_ij 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 paw_ij treated by current proc if not present, will be calculated in the present routine mpi_atmtab_out= --optional-- indexes of the output paw_ij 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
[paw_ij_out(:)]<type(paw_ij_type)>= --optional-- if present, the redistributed datastructure does not replace the input one but is delivered in paw_ij_out if not present, input and output datastructure are the same.
SIDE EFFECTS
paw_ij(:)<type(paw_ij_type)>= input (and eventually output) paw_ij datastructures
PARENTS
m_paral_pert
CHILDREN
SOURCE
2079 subroutine paw_ij_redistribute(paw_ij,mpi_comm_in,mpi_comm_out,& 2080 & natom,mpi_atmtab_in,mpi_atmtab_out,paw_ij_out,& 2081 & SendAtomProc,SendAtomList,RecvAtomProc,RecvAtomList) 2082 2083 2084 !This section has been created automatically by the script Abilint (TD). 2085 !Do not modify the following lines by hand. 2086 #undef ABI_FUNC 2087 #define ABI_FUNC 'paw_ij_redistribute' 2088 !End of the abilint section 2089 2090 implicit none 2091 2092 !Arguments ------------------------------------ 2093 !scalars 2094 integer,intent(in) :: mpi_comm_in,mpi_comm_out 2095 integer,optional,intent(in) :: natom 2096 !arrays 2097 integer,intent(in),optional,target :: mpi_atmtab_in(:),mpi_atmtab_out(:) 2098 type(paw_ij_type),allocatable,intent(inout) :: paw_ij(:) 2099 type(paw_ij_type),pointer,optional :: paw_ij_out(:) !vz_i 2100 integer,intent(in),optional :: SendAtomProc(:),SendAtomList(:),RecvAtomProc(:),RecvAtomList(:) 2101 2102 !Local variables------------------------------- 2103 !scalars 2104 integer :: algo_option,i1,iatom,iat_in,iat_out,ierr,iircv,iisend,imsg,imsg_current 2105 integer :: imsg1,iproc_rcv,iproc_send,ireq,me_exch,mpi_comm_exch,my_natom_in,my_natom_out,my_tag,natom_tot 2106 integer :: nb_dp,nb_int,nb_msg,nbmsg_incoming,nbrecv,nbrecvmsg,nbsendreq,nbsent,nbsend,next,npaw_ij_sent 2107 integer :: nproc_in,nproc_out 2108 logical :: flag,in_place,message_yet_prepared,my_atmtab_in_allocated,my_atmtab_out_allocated,paral_atom 2109 !arrays 2110 integer :: buf_size(3),request1(3) 2111 integer,pointer :: my_atmtab_in(:),my_atmtab_out(:) 2112 integer,allocatable :: atmtab_send(:),atm_indx_in(:),atm_indx_out(:),From(:),buf_int1(:),request(:) 2113 integer,allocatable,target:: buf_int(:) 2114 integer,pointer :: buf_ints(:) 2115 logical,allocatable :: msg_pick(:) 2116 real(dp),allocatable :: buf_dp1(:) 2117 real(dp),allocatable,target :: buf_dp(:) 2118 real(dp),pointer :: buf_dps(:) 2119 type(coeffi1_type),target,allocatable :: tab_buf_int(:),tab_buf_atom(:) 2120 type(coeff1_type),target,allocatable :: tab_buf_dp(:) 2121 type(paw_ij_type),allocatable :: paw_ij_all(:) 2122 type(paw_ij_type),pointer :: paw_ij_out1(:) 2123 2124 ! ************************************************************************* 2125 2126 !@Paw_ij_type 2127 2128 in_place=(.not.present(paw_ij_out)) 2129 my_natom_in=size(paw_ij) 2130 2131 !If not "in_place", destroy the output datastructure 2132 if (.not.in_place) then 2133 if (associated(paw_ij_out)) then 2134 call paw_ij_free(paw_ij_out) 2135 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_out) 2136 end if 2137 end if 2138 2139 !Special sequential case 2140 if (mpi_comm_in==xmpi_comm_self.and.mpi_comm_out==xmpi_comm_self) then 2141 if ((.not.in_place).and.(my_natom_in>0)) then 2142 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_in)) 2143 call paw_ij_nullify(paw_ij_out) 2144 call paw_ij_copy(paw_ij,paw_ij_out) 2145 end if 2146 return 2147 end if 2148 2149 !Get total natom 2150 if (present(natom)) then 2151 natom_tot=natom 2152 else 2153 natom_tot=my_natom_in 2154 call xmpi_sum(natom_tot,mpi_comm_in,ierr) 2155 end if 2156 2157 !Select input distribution 2158 if (present(mpi_atmtab_in)) then 2159 my_atmtab_in => mpi_atmtab_in 2160 my_atmtab_in_allocated=.false. 2161 else 2162 call get_my_atmtab(mpi_comm_in,my_atmtab_in,my_atmtab_in_allocated,& 2163 & paral_atom,natom_tot,my_natom_in) 2164 end if 2165 2166 !Select output distribution 2167 if (present(mpi_atmtab_out)) then 2168 my_natom_out=size(mpi_atmtab_out) 2169 my_atmtab_out => mpi_atmtab_out 2170 my_atmtab_out_allocated=.false. 2171 else 2172 call get_my_natom(mpi_comm_out,my_natom_out,natom_tot) 2173 call get_my_atmtab(mpi_comm_out,my_atmtab_out,my_atmtab_out_allocated,& 2174 & paral_atom,natom_tot) 2175 end if 2176 2177 !Select algo according to optional input arguments 2178 algo_option=1 2179 if (present(SendAtomProc).and.present(SendAtomList).and.& 2180 & present(RecvAtomProc).and.present(RecvAtomList)) algo_option=2 2181 2182 2183 !Brute force algorithm (allgather + scatter) 2184 !--------------------------------------------------------- 2185 if (algo_option==1) then 2186 2187 LIBPAW_DATATYPE_ALLOCATE(paw_ij_all,(natom_tot)) 2188 call paw_ij_nullify(paw_ij_all) 2189 call paw_ij_copy(paw_ij,paw_ij_all,comm_atom=mpi_comm_in,mpi_atmtab=my_atmtab_in) 2190 if (in_place) then 2191 call paw_ij_free(paw_ij) 2192 LIBPAW_DATATYPE_DEALLOCATE(paw_ij) 2193 LIBPAW_DATATYPE_ALLOCATE(paw_ij,(my_natom_out)) 2194 call paw_ij_nullify(paw_ij) 2195 call paw_ij_copy(paw_ij_all,paw_ij,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 2196 else 2197 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_out)) 2198 call paw_ij_nullify(paw_ij_out) 2199 call paw_ij_copy(paw_ij_all,paw_ij_out,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 2200 end if 2201 call paw_ij_free(paw_ij_all) 2202 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_all) 2203 2204 2205 !Asynchronous algorithm (asynchronous communications) 2206 !--------------------------------------------------------- 2207 else if (algo_option==2) then 2208 2209 nbsend=size(SendAtomProc) ; nbrecv=size(RecvAtomProc) 2210 2211 if (in_place) then 2212 if (my_natom_out > 0) then 2213 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out1,(my_natom_out)) 2214 call paw_ij_nullify(paw_ij_out1) 2215 else 2216 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out1,(0)) 2217 end if 2218 else 2219 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_out)) 2220 call paw_ij_nullify(paw_ij_out) 2221 paw_ij_out1=>paw_ij_out 2222 end if 2223 2224 nproc_in=xmpi_comm_size(mpi_comm_in) 2225 nproc_out=xmpi_comm_size(mpi_comm_out) 2226 if (nproc_in<=nproc_out) mpi_comm_exch=mpi_comm_out 2227 if (nproc_in>nproc_out) mpi_comm_exch=mpi_comm_in 2228 me_exch=xmpi_comm_rank(mpi_comm_exch) 2229 2230 2231 ! Dimension put to the maximum to send 2232 LIBPAW_ALLOCATE(atmtab_send,(nbsend)) 2233 LIBPAW_ALLOCATE(atm_indx_in,(natom_tot)) 2234 atm_indx_in=-1 2235 do iatom=1,my_natom_in 2236 atm_indx_in(my_atmtab_in(iatom))=iatom 2237 end do 2238 LIBPAW_ALLOCATE(atm_indx_out,(natom_tot)) 2239 atm_indx_out=-1 2240 do iatom=1,my_natom_out 2241 atm_indx_out(my_atmtab_out(iatom))=iatom 2242 end do 2243 2244 LIBPAW_DATATYPE_ALLOCATE(tab_buf_int,(nbsend)) 2245 LIBPAW_DATATYPE_ALLOCATE(tab_buf_dp,(nbsend)) 2246 LIBPAW_DATATYPE_ALLOCATE(tab_buf_atom,(nbsend)) 2247 LIBPAW_ALLOCATE(request,(3*nbsend)) 2248 2249 ! A send buffer in an asynchrone communication couldn't be deallocate before it has been receive 2250 nbsent=0 ; ireq=0 ; iisend=0 ; nbsendreq=0 ; nb_msg=0 2251 do iisend=1,nbsend 2252 iproc_rcv=SendAtomProc(iisend) 2253 next=-1 2254 if (iisend < nbsend) next=SendAtomProc(iisend+1) 2255 if (iproc_rcv /= me_exch) then 2256 nbsent=nbsent+1 2257 atmtab_send(nbsent)=SendAtomList(iisend) ! we groups the atoms sends to the same process 2258 if (iproc_rcv /= next) then 2259 if (nbsent > 0 ) then 2260 ! Check if message has been yet prepared 2261 message_yet_prepared=.false. 2262 do imsg=1,nb_msg 2263 if (size(tab_buf_atom(imsg)%value) /= nbsent) then 2264 cycle 2265 else 2266 do imsg1=1,nbsent 2267 if (tab_buf_atom(imsg)%value(imsg1)/=atmtab_send(imsg1)) exit 2268 message_yet_prepared=.true. 2269 imsg_current=imsg 2270 end do 2271 end if 2272 enddo 2273 ! Create the message 2274 if (.not.message_yet_prepared) then 2275 nb_msg=nb_msg+1 2276 call paw_ij_isendreceive_fillbuffer( & 2277 & paw_ij,atmtab_send,atm_indx_in,nbsent,buf_int,nb_int,buf_dp,nb_dp) 2278 LIBPAW_ALLOCATE(tab_buf_int(nb_msg)%value,(nb_int)) 2279 LIBPAW_ALLOCATE(tab_buf_dp(nb_msg)%value,(nb_dp)) 2280 tab_buf_int(nb_msg)%value(1:nb_int)=buf_int(1:nb_int) 2281 tab_buf_dp(nb_msg)%value(1:nb_dp)=buf_dp(1:nb_dp) 2282 LIBPAW_DEALLOCATE(buf_int) 2283 LIBPAW_DEALLOCATE(buf_dp) 2284 LIBPAW_ALLOCATE(tab_buf_atom(nb_msg)%value, (nbsent)) 2285 tab_buf_atom(nb_msg)%value(1:nbsent)=atmtab_send(1:nbsent) 2286 imsg_current=nb_msg 2287 end if 2288 ! Communicate 2289 buf_size(1)=size(tab_buf_int(imsg_current)%value) 2290 buf_size(2)=size(tab_buf_dp(imsg_current)%value) 2291 buf_size(3)=nbsent 2292 buf_ints=>tab_buf_int(imsg_current)%value 2293 buf_dps=>tab_buf_dp(imsg_current)%value 2294 my_tag=200 2295 ireq=ireq+1 2296 call xmpi_isend(buf_size,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2297 my_tag=201 2298 ireq=ireq+1 2299 call xmpi_isend(buf_ints,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2300 my_tag=202 2301 ireq=ireq+1 2302 call xmpi_isend(buf_dps,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2303 nbsendreq=ireq 2304 nbsent=0 2305 end if 2306 end if 2307 else ! Just a renumbering, not a sending 2308 iat_in=atm_indx_in(SendAtomList(iisend)) 2309 iat_out=atm_indx_out(my_atmtab_in(iat_in)) 2310 call paw_ij_copy(paw_ij(iat_in:iat_in),paw_ij_out1(iat_out:iat_out)) 2311 nbsent=0 2312 end if 2313 end do 2314 2315 LIBPAW_ALLOCATE(From,(nbrecv)) 2316 From(:)=-1 ; nbrecvmsg=0 2317 do iircv=1,nbrecv 2318 iproc_send=RecvAtomProc(iircv) !receive from (RcvAtomProc is sorted by growing process) 2319 next=-1 2320 if (iircv < nbrecv) next=RecvAtomProc(iircv+1) 2321 if (iproc_send /= me_exch .and. iproc_send/=next) then 2322 nbrecvmsg=nbrecvmsg+1 2323 From(nbrecvmsg)=iproc_send 2324 end if 2325 end do 2326 2327 LIBPAW_ALLOCATE(msg_pick,(nbrecvmsg)) 2328 msg_pick=.false. 2329 nbmsg_incoming=nbrecvmsg 2330 do while (nbmsg_incoming > 0) 2331 do i1=1,nbrecvmsg 2332 if (.not.msg_pick(i1)) then 2333 iproc_send=From(i1) 2334 flag=.false. 2335 my_tag=200 2336 call xmpi_iprobe(iproc_send,my_tag,mpi_comm_exch,flag,ierr) 2337 if (flag) then 2338 msg_pick(i1)=.true. 2339 call xmpi_irecv(buf_size,iproc_send,my_tag,mpi_comm_exch,request1(1),ierr) 2340 call xmpi_wait(request1(1),ierr) 2341 nb_int=buf_size(1) 2342 nb_dp=buf_size(2) 2343 npaw_ij_sent=buf_size(3) 2344 LIBPAW_ALLOCATE(buf_int1,(nb_int)) 2345 LIBPAW_ALLOCATE(buf_dp1,(nb_dp)) 2346 my_tag=201 2347 call xmpi_irecv(buf_int1,iproc_send,my_tag,mpi_comm_exch,request1(2),ierr) 2348 my_tag=202 2349 call xmpi_irecv(buf_dp1,iproc_send,my_tag,mpi_comm_exch,request1(3),ierr) 2350 call xmpi_waitall(request1(2:3),ierr) 2351 call paw_ij_isendreceive_getbuffer(paw_ij_out1,npaw_ij_sent,atm_indx_out,buf_int1,buf_dp1) 2352 nbmsg_incoming=nbmsg_incoming-1 2353 LIBPAW_DEALLOCATE(buf_int1) 2354 LIBPAW_DEALLOCATE(buf_dp1) 2355 end if 2356 end if 2357 end do 2358 end do 2359 LIBPAW_DEALLOCATE(msg_pick) 2360 2361 if (in_place) then 2362 call paw_ij_free(paw_ij) 2363 LIBPAW_DATATYPE_DEALLOCATE(paw_ij) 2364 LIBPAW_DATATYPE_ALLOCATE(paw_ij,(my_natom_out)) 2365 call paw_ij_copy(paw_ij_out1,paw_ij) 2366 call paw_ij_free(paw_ij_out1) 2367 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_out1) 2368 end if 2369 2370 ! Wait for deallocating arrays that all sending operations has been realized 2371 if (nbsendreq > 0) then 2372 call xmpi_waitall(request(1:nbsendreq),ierr) 2373 end if 2374 2375 ! Deallocate buffers 2376 do i1=1,nb_msg 2377 LIBPAW_DEALLOCATE(tab_buf_int(i1)%value) 2378 LIBPAW_DEALLOCATE(tab_buf_dp(i1)%value) 2379 LIBPAW_DEALLOCATE(tab_buf_atom(i1)%value) 2380 end do 2381 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_int) 2382 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_dp) 2383 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_atom) 2384 LIBPAW_DEALLOCATE(From) 2385 LIBPAW_DEALLOCATE(request) 2386 LIBPAW_DEALLOCATE(atmtab_send) 2387 LIBPAW_DEALLOCATE(atm_indx_in) 2388 LIBPAW_DEALLOCATE(atm_indx_out) 2389 2390 end if !algo_option 2391 2392 !Eventually release temporary pointers 2393 call free_my_atmtab(my_atmtab_in,my_atmtab_in_allocated) 2394 call free_my_atmtab(my_atmtab_out,my_atmtab_out_allocated) 2395 2396 end subroutine paw_ij_redistribute
m_paw_ij/paw_ij_reset_flags [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
paw_ij_reset_flags
FUNCTION
Set all paw_ij flags to 1 (force the recomputation of all arrays)
SIDE EFFECTS
Paw_ij<type(paw_ij_type)>=paw_ij structure
PARENTS
d2frnl,dfpt_nstpaw,dfpt_scfcv,scfcv
CHILDREN
SOURCE
2418 subroutine paw_ij_reset_flags(Paw_ij,all,dijhartree,self_consistent) 2419 2420 2421 !This section has been created automatically by the script Abilint (TD). 2422 !Do not modify the following lines by hand. 2423 #undef ABI_FUNC 2424 #define ABI_FUNC 'paw_ij_reset_flags' 2425 !End of the abilint section 2426 2427 implicit none 2428 2429 !Arguments ------------------------------------ 2430 !scalars 2431 logical,optional,intent(in) :: all,dijhartree,self_consistent 2432 !arrays 2433 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 2434 2435 !Local variables------------------------------- 2436 integer :: iat,natom 2437 logical :: all_,dijhartree_,self_consistent_ 2438 2439 ! ************************************************************************* 2440 2441 !@Paw_ij_type 2442 2443 natom=SIZE(Paw_ij);if (natom==0) return 2444 all_=.true.;if (present(all)) all_=all 2445 dijhartree_=.false.;if (present(dijhartree)) dijhartree_=dijhartree 2446 self_consistent_=.false.;if (present(self_consistent)) self_consistent_=self_consistent 2447 2448 if (dijhartree_) then 2449 do iat=1,natom 2450 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2451 end do 2452 2453 else if (self_consistent_) then 2454 do iat=1,natom 2455 if (Paw_ij(iat)%has_dij >0) Paw_ij(iat)%has_dij =1 2456 if (Paw_ij(iat)%has_dijexxc >0) Paw_ij(iat)%has_dijexxc =1 2457 if (Paw_ij(iat)%has_dijfock >0) Paw_ij(iat)%has_dijfock =1 2458 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2459 if (Paw_ij(iat)%has_dijhat >0) Paw_ij(iat)%has_dijhat =1 2460 if (Paw_ij(iat)%has_dijnd >0) Paw_ij(iat)%has_dijnd =1 2461 if (Paw_ij(iat)%has_dijso >0) Paw_ij(iat)%has_dijso =1 2462 if (Paw_ij(iat)%has_dijU >0) Paw_ij(iat)%has_dijU =1 2463 if (Paw_ij(iat)%has_dijxc >0) Paw_ij(iat)%has_dijxc =1 2464 if (Paw_ij(iat)%has_dijxc_hat >0) Paw_ij(iat)%has_dijxc_hat =1 2465 if (Paw_ij(iat)%has_dijxc_val >0) Paw_ij(iat)%has_dijxc_val =1 2466 if (Paw_ij(iat)%has_exexch_pot>0) Paw_ij(iat)%has_exexch_pot=1 2467 if (Paw_ij(iat)%has_pawu_occ >0) Paw_ij(iat)%has_pawu_occ =1 2468 end do 2469 2470 else if (all_) then 2471 do iat=1,natom 2472 if (Paw_ij(iat)%has_dij >0) Paw_ij(iat)%has_dij =1 2473 if (Paw_ij(iat)%has_dij0 >0) Paw_ij(iat)%has_dij0 =1 2474 if (Paw_ij(iat)%has_dijexxc >0) Paw_ij(iat)%has_dijexxc =1 2475 if (Paw_ij(iat)%has_dijfock >0) Paw_ij(iat)%has_dijfock =1 2476 if (Paw_ij(iat)%has_dijfr >0) Paw_ij(iat)%has_dijfr =1 2477 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2478 if (Paw_ij(iat)%has_dijhat >0) Paw_ij(iat)%has_dijhat =1 2479 if (Paw_ij(iat)%has_dijnd >0) Paw_ij(iat)%has_dijnd =1 2480 if (Paw_ij(iat)%has_dijso >0) Paw_ij(iat)%has_dijso =1 2481 if (Paw_ij(iat)%has_dijU >0) Paw_ij(iat)%has_dijU =1 2482 if (Paw_ij(iat)%has_dijxc >0) Paw_ij(iat)%has_dijxc =1 2483 if (Paw_ij(iat)%has_dijxc_hat >0) Paw_ij(iat)%has_dijxc_hat =1 2484 if (Paw_ij(iat)%has_dijxc_val >0) Paw_ij(iat)%has_dijxc_val =1 2485 if (Paw_ij(iat)%has_exexch_pot>0) Paw_ij(iat)%has_exexch_pot=1 2486 if (Paw_ij(iat)%has_pawu_occ >0) Paw_ij(iat)%has_pawu_occ =1 2487 end do 2488 end if 2489 2490 end subroutine paw_ij_reset_flags
m_paw_ij/paw_ij_type [ Types ]
[ Top ] [ m_paw_ij ] [ Types ]
NAME
paw_ij_type
FUNCTION
For PAW, various arrays given on (i,j) (partial waves) channels
SOURCE
53 type,public :: paw_ij_type 54 55 !Integer scalars 56 57 integer :: cplex_rf 58 ! cplex_rf=1 if on-site PAW quantities are real 59 ! cplex_rf=2 if on-site PAW quantities are complex due to response at q<>0 (epx(-iqr) phase) 60 61 integer :: cplex_dij 62 ! cplex_dij=1 if dij are real 63 ! cplex_dij=2 if dij are complex (spin-orbit, non-collinear magnetism, magnetic field, ...) 64 65 integer :: has_dij=0 66 ! 1 if dij is allocated 67 ! 2 if dij is already computed 68 69 integer :: has_dij0=0 70 ! 1 if dij0 is allocated 71 ! 2 if dij0 is already computed 72 73 integer :: has_dijexxc=0 74 ! 1 if dijexxc is associated and used, 0 otherwise 75 ! 2 if dijexxc is already computed 76 77 integer :: has_dijfock=0 78 ! 1 if dijfock is allocated 79 ! 2 if dijfock is already computed 80 81 integer :: has_dijfr=0 82 ! 1 if dijfr is allocated 83 ! 2 if dijfr is already computed 84 85 integer :: has_dijhartree=0 86 ! 1 if dijhartree is allocated 87 ! 2 if dijhartree is already computed 88 89 integer :: has_dijhat=0 90 ! 1 if dijhat is allocated 91 ! 2 if dijhat is already computed 92 93 integer :: has_dijnd=0 94 ! on site term due to nuclear dipole moment 95 ! 1 if dijnd is associated and used, 0 otherwise 96 ! 2 if dijnd is already computed 97 98 integer :: has_dijso=0 99 ! 1 if dijso is associated and used, 0 otherwise 100 ! 2 if dijso is already computed 101 102 integer :: has_dijU=0 103 ! 1 if dijU is associated and used, 0 otherwise 104 ! 2 if dijU is already computed 105 106 integer :: has_dijxc=0 107 ! 1 if dijxc is associated and used, 0 otherwise 108 ! 2 if dijxc is already computed 109 110 integer :: has_dijxc_hat=0 111 ! 1 if dijxc_hat is associated and used, 0 otherwise 112 ! 2 if dijxc_hat is already computed 113 114 integer :: has_dijxc_val=0 115 ! 1 if dijxc_val is associated and used, 0 otherwise 116 ! 2 if dijxc_val is already computed 117 118 integer :: has_exexch_pot=0 119 ! 1 if PAW+(local exact exchange) potential is allocated 120 121 integer :: has_pawu_occ=0 122 ! 1 if PAW+U occupations are allocated 123 124 integer :: itypat 125 ! itypat=type of the atom 126 127 integer :: lmn_size 128 ! Number of (l,m,n) elements for the paw basis 129 130 integer :: lmn2_size 131 ! lmn2_size=lmn_size*(lmn_size+1)/2 132 ! where lmn_size is the number of (l,m,n) elements for the paw basis 133 134 integer :: ndij 135 ! Number of components of dij 136 ! Usually ndij=nspden, except for nspinor==2 (where ndij=nspinor**2) 137 138 integer :: nspden 139 ! Number of spin-density components (may be different from dtset%nspden if spin-orbit) 140 141 integer :: nsppol 142 ! Number of independant spin-components 143 144 145 !Real (real(dp)) arrays 146 147 real(dp), allocatable :: dij(:,:) 148 ! dij(cplex_rf*cplex_dij*lmn2_size,ndij) 149 ! Dij term (non-local operator) 150 ! May be complex if cplex_dij=2 or cplex_rf=2 151 ! Storage for the 1st dimension: 152 ! For each klmn=(i,j), (real,imaginary) parts 153 ! If cplex_rf=2, first all values for iplex_rf=1 (real part of exp(-iqr)), 154 ! then all values for iplex_rf=2 (imaginary part if exp(-iqr)) 155 ! Storage for the 2st dimension: 156 ! dij(:,1) contains Dij^up-up 157 ! dij(:,2) contains Dij^dn-dn 158 ! dij(:,3) contains Dij^up-dn (only if nspinor=2) 159 ! dij(:,4) contains Dij^dn-up (only if nspinor=2) 160 161 real(dp), allocatable :: dij0(:) 162 ! dij0(lmn2_size) 163 ! Atomic part of Dij (read from PAW dataset) 164 ! Same storage as Dij (see above); always real, spin independent 165 166 real(dp), allocatable :: dijexxc(:,:) 167 ! dijexxc(cplex_dij*lmn2_size,ndij) 168 ! On-site matrix elements of the Fock operator (Local Exact exchange implementation) 169 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 170 171 real(dp), allocatable :: dijfock(:,:) 172 ! dijfock(cplex_dij*lmn2_size,ndij) 173 ! Dij_fock term 174 ! Contains all contributions to Dij from Fock exchange 175 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 176 177 real(dp), allocatable :: dijfr(:,:) 178 ! dijhat(cplex_rf*cplex_dij*lmn2_size,ndij) 179 ! For response function calculation only 180 ! RF Frozen part of Dij (depends on q vector but not on 1st-order wave function) 181 ! Same storage as Dij (see above) 182 183 real(dp), allocatable :: dijhartree(:) 184 ! dijhartree(cplex_rf*lmn2_size) 185 ! Dij_hartree term; contains all contributions to Dij from hartree 186 ! Warning: Dimensioned only by cplex_rf (exp(-iqr)), not cplex_dij 187 ! Same storage as Dij (see above); spin independent 188 189 real(dp), allocatable :: dijhat(:,:) 190 ! dijhat(cplex_rf*cplex_dij*lmn2_size,ndij) 191 ! Dij_hat term (non-local operator) i.e \sum_LM \int_FFT Q_{ij}^{LM} vtrial 192 ! Same storage as Dij (see above) 193 ! Same storage as Dij (see above) 194 195 real(dp), allocatable :: dijnd(:,:) 196 ! dijnd(cplex_dij*lmn2_size,ndij) 197 ! On-site matrix elements of -\frac{1}{c}\mu\cdot L/r^3 198 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 199 200 real(dp), allocatable :: dijso(:,:) 201 ! dijso(cplex_rf*cplex_dij*lmn2_size,ndij) 202 ! On-site matrix elements of L.S i.e <phi_i|L.S|phi_j> 203 ! Same storage as Dij (see above) 204 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 205 206 real(dp), allocatable :: dijU(:,:) 207 ! dijU(cplex_rf*cplex_dij*lmn2_size,ndij) 208 ! On-site matrix elements of the U part of the PAW Hamiltonian. 209 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 210 211 real(dp), allocatable :: dijxc(:,:) 212 ! dijxc(cplex_rf*cplex_dij*lmn2_size,ndij) 213 ! On-site matrix elements of vxc i.e 214 ! <phi_i|vxc[n1+nc]|phi_j> - <tphi_i|vxc(tn1+nhat+tnc]|tphi_j> 215 ! Same storage as Dij (see above) 216 217 real(dp), allocatable :: dijxc_hat(:,:) 218 ! dijxc_hat(cplex_dij*lmn2_size,ndij) 219 ! Dij_hat term i.e \sum_LM \int_FFT Q_{ij}^{LM} Vxc 220 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 221 222 real(dp), allocatable :: dijxc_val(:,:) 223 ! dijxc_val(cplex_dij*lmn2_size,ndij) 224 ! Onsite matrix elements of valence-only vxc i.e 225 ! <phi_i|vxc[n1]|phi_j> - <tphi_i|vxc(tn1+nhat]|tphi_j> 226 ! Same storage as Dij (see above); not available for RF (cplex_rf=1) 227 228 real(dp), allocatable :: noccmmp(:,:,:,:) 229 ! noccmmp(cplex_dij,2*lpawu+1,2*lpawu+1,nocc_nspden) 230 ! cplex_dij=1 if collinear 231 ! cplex_dij=2 if spin orbit is used 232 ! cplex_dij=2 is used if non-collinear (for coherence, it is not necessary in this case, however) 233 ! gives occupation matrix for lda+u (computed in setnoccmmp) 234 ! Stored as: noccmmp(:,:,1)= n^{up,up}_{m,mp} 235 ! noccmmp(:,:,2)= n^{dn,dn}_{m,mp} 236 ! noccmmp(:,:,3)= n^{up,dn}_{m,mp} 237 ! noccmmp(:,:,4)= n^{dn,up}_{m,mp} 238 ! noccmmp(m,mp,:) is computed from rhoij(klmn) with m=klmntomn(2)>mp=klmntomn(1) 239 240 real(dp), allocatable :: nocctot(:) 241 ! nocctot(ndij) 242 ! gives trace of occupation matrix for lda+u (computed in pawdenpot) 243 ! for each value of ispden (1 or 2) 244 245 real(dp), allocatable :: vpawx(:,:,:) 246 ! vpawx(2*lexexch+1,2*lexexch+1,nspden) 247 ! exact exchange potential 248 249 end type paw_ij_type 250 251 !public procedures. 252 public :: paw_ij_init ! Creation method 253 public :: paw_ij_free ! Free memory 254 public :: paw_ij_nullify 255 public :: paw_ij_copy ! Copy object 256 public :: paw_ij_print ! Printout of the object 257 public :: paw_ij_gather ! MPI gather 258 public :: paw_ij_redistribute ! MPI redistribute 259 public :: paw_ij_reset_flags ! Resent the internal flags. 260 261 !private procedures. 262 private :: paw_ij_isendreceive_getbuffer 263 private :: paw_ij_isendreceive_fillbuffer