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-2024 ABINIT group (MT, FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
23 #include "libpaw.h" 24 25 MODULE m_paw_ij 26 27 USE_DEFS 28 USE_MSG_HANDLING 29 USE_MPI_WRAPPERS 30 USE_MEMORY_PROFILING 31 32 use m_paral_atom, only : get_my_atmtab, free_my_atmtab, get_my_natom 33 use m_pawtab, only : pawtab_type 34 use m_paw_io, only : pawio_print_ij 35 36 implicit none 37 38 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
SOURCE
2641 subroutine paw_ij_isendreceive_fillbuffer(paw_ij,atmtab_send,atm_indx_send,npaw_ij_send,& 2642 & buf_int,buf_int_size,buf_dp,buf_dp_size) 2643 2644 !Arguments ------------------------------------ 2645 !scalars 2646 integer,intent(out) :: buf_int_size,buf_dp_size 2647 integer,intent(in) :: npaw_ij_send 2648 !arrays 2649 integer,intent(in) :: atmtab_send(:),atm_indx_send(:) 2650 integer,allocatable,intent(out) :: buf_int(:) 2651 real(dp),allocatable,intent(out):: buf_dp(:) 2652 type(paw_ij_type),target,intent(in) :: paw_ij(:) 2653 2654 !Local variables------------------------------- 2655 !scalars 2656 integer :: cplxdij_lmn2_size,cplxdijq_lmn2_size,cplxq_lmn2_size 2657 integer :: iatom_tot,ii,ij,indx_dp,indx_int,ipaw_ij_send 2658 integer :: lmn2_size,ndij,nocc,nspden,sz1,sz2,sz3 2659 character(len=500) :: msg 2660 type(Paw_ij_type),pointer :: paw_ij1 2661 !arrays 2662 2663 ! ********************************************************************* 2664 2665 !Compute sizes of buffers 2666 buf_int_size=0;buf_dp_size=0 2667 do ipaw_ij_send=1,npaw_ij_send 2668 iatom_tot=atmtab_send(ipaw_ij_send) 2669 ij = atm_indx_send(iatom_tot) 2670 paw_ij1=>paw_ij(ij) 2671 lmn2_size=paw_ij1%lmn2_size 2672 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2673 cplxq_lmn2_size=paw_ij1%qphase*lmn2_size 2674 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij1%qphase 2675 ndij=paw_ij1%ndij 2676 buf_int_size=buf_int_size+24 2677 if (paw_ij1%has_dij==2) then 2678 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2679 end if 2680 if (paw_ij1%has_dij0==2) then 2681 buf_dp_size=buf_dp_size +lmn2_size 2682 end if 2683 if (paw_ij1%has_dijexxc==2) then 2684 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2685 end if 2686 if (paw_ij1%has_dijfock==2) then 2687 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2688 end if 2689 if (paw_ij1%has_dijfr==2) then 2690 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2691 end if 2692 if (paw_ij1%has_dijhartree==2) then 2693 buf_dp_size=buf_dp_size +cplxq_lmn2_size 2694 end if 2695 if (paw_ij1%has_dijhat==2) then 2696 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2697 end if 2698 if (paw_ij1%has_dijnd==2) then 2699 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2700 end if 2701 if (paw_ij1%has_dijso==2) then 2702 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2703 end if 2704 if (paw_ij1%has_dijU==2) then 2705 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2706 end if 2707 if (paw_ij1%has_dijxc==2) then 2708 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 2709 end if 2710 if (paw_ij1%has_dijxc_hat==2) then 2711 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2712 end if 2713 if (paw_ij1%has_dijxc_val==2) then 2714 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 2715 end if 2716 if (paw_ij1%has_pawu_occ>=1) then 2717 buf_int_size=buf_int_size+4 2718 buf_dp_size=buf_dp_size & 2719 & +size(paw_ij1%nocctot) & 2720 & +size(paw_ij1%noccmmp) 2721 end if 2722 if (paw_ij1%has_exexch_pot>=1) then 2723 buf_int_size=buf_int_size+3 2724 buf_dp_size=buf_dp_size +size(paw_ij1%vpawx) 2725 end if 2726 end do 2727 2728 !Fill input buffers 2729 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 2730 LIBPAW_ALLOCATE(buf_dp,(buf_dp_size)) 2731 indx_int=1;indx_dp=1 2732 do ipaw_ij_send=1,npaw_ij_send 2733 iatom_tot=atmtab_send(ipaw_ij_send) 2734 ij = atm_indx_send(iatom_tot) 2735 paw_ij1=>paw_ij(ij) 2736 nspden=paw_ij1%nspden 2737 buf_int(indx_int)=iatom_tot ;indx_int=indx_int+1 2738 buf_int(indx_int)=paw_ij1%cplex_dij ;indx_int=indx_int+1 2739 buf_int(indx_int)=paw_ij1%qphase ;indx_int=indx_int+1 2740 buf_int(indx_int)=paw_ij1%itypat ;indx_int=indx_int+1 2741 buf_int(indx_int)=nspden ;indx_int=indx_int+1 2742 buf_int(indx_int)=paw_ij1%nsppol ;indx_int=indx_int+1 2743 buf_int(indx_int)=paw_ij1%lmn_size ;indx_int=indx_int+1 2744 buf_int(indx_int)=paw_ij1%lmn2_size ;indx_int=indx_int+1 2745 buf_int(indx_int)=paw_ij1%ndij ;indx_int=indx_int+1 2746 buf_int(indx_int)=paw_ij1%has_dij ;indx_int=indx_int+1 2747 buf_int(indx_int)=paw_ij1%has_dij0 ;indx_int=indx_int+1 2748 buf_int(indx_int)=paw_ij1%has_dijexxc ;indx_int=indx_int+1 2749 buf_int(indx_int)=paw_ij1%has_dijfock ;indx_int=indx_int+1 2750 buf_int(indx_int)=paw_ij1%has_dijfr ;indx_int=indx_int+1 2751 buf_int(indx_int)=paw_ij1%has_dijhartree ;indx_int=indx_int+1 2752 buf_int(indx_int)=paw_ij1%has_dijhat ;indx_int=indx_int+1 2753 buf_int(indx_int)=paw_ij1%has_dijnd ;indx_int=indx_int+1 2754 buf_int(indx_int)=paw_ij1%has_dijso ;indx_int=indx_int+1 2755 buf_int(indx_int)=paw_ij1%has_dijU ;indx_int=indx_int+1 2756 buf_int(indx_int)=paw_ij1%has_dijxc ;indx_int=indx_int+1 2757 buf_int(indx_int)=paw_ij1%has_dijxc_hat ;indx_int=indx_int+1 2758 buf_int(indx_int)=paw_ij1%has_dijxc_val ;indx_int=indx_int+1 2759 buf_int(indx_int)=paw_ij1%has_exexch_pot ;indx_int=indx_int+1 2760 buf_int(indx_int)=paw_ij1%has_pawu_occ ;indx_int=indx_int+1 2761 lmn2_size=paw_ij1%lmn2_size 2762 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2763 cplxq_lmn2_size=paw_ij1%qphase*lmn2_size 2764 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij1%qphase 2765 ndij=paw_ij1%ndij 2766 if (paw_ij1%has_dij==2) then 2767 ii=cplxdijq_lmn2_size*ndij 2768 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dij,(/ii/)) 2769 indx_dp=indx_dp+ii 2770 end if 2771 if (paw_ij1%has_dij0==2) then 2772 ii=lmn2_size 2773 buf_dp(indx_dp:indx_dp+lmn2_size-1)=paw_ij1%dij0(:) 2774 indx_dp=indx_dp+lmn2_size 2775 end if 2776 if (paw_ij1%has_dijexxc==2) then 2777 ii=cplxdij_lmn2_size*ndij 2778 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijexxc,(/ii/)) 2779 indx_dp=indx_dp+ii 2780 end if 2781 if (paw_ij1%has_dijfock==2) then 2782 ii=cplxdij_lmn2_size*ndij 2783 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijfock,(/ii/)) 2784 indx_dp=indx_dp+ii 2785 end if 2786 if (paw_ij1%has_dijfr==2) then 2787 ii=cplxdijq_lmn2_size*ndij 2788 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijfr,(/ii/)) 2789 indx_dp=indx_dp+ii 2790 end if 2791 if (paw_ij1%has_dijhartree==2) then 2792 ii=cplxq_lmn2_size 2793 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij1%dijhartree(:) 2794 indx_dp=indx_dp+ii 2795 end if 2796 if (paw_ij1%has_dijhat==2) then 2797 ii=cplxdijq_lmn2_size*ndij 2798 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijhat,(/ii/)) 2799 indx_dp=indx_dp+ii 2800 end if 2801 if (paw_ij1%has_dijnd==2) then 2802 ii=cplxdij_lmn2_size*ndij 2803 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijnd,(/ii/)) 2804 indx_dp=indx_dp+ii 2805 end if 2806 if (paw_ij1%has_dijso==2) then 2807 ii=cplxdijq_lmn2_size*ndij 2808 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijso,(/ii/)) 2809 indx_dp=indx_dp+ii 2810 end if 2811 if (paw_ij1%has_dijU==2) then 2812 ii=cplxdijq_lmn2_size*ndij 2813 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijU,(/ii/)) 2814 indx_dp=indx_dp+ii 2815 end if 2816 if (paw_ij1%has_dijxc==2) then 2817 ii=cplxdijq_lmn2_size*ndij 2818 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc,(/ii/)) 2819 indx_dp=indx_dp+ii 2820 end if 2821 if (paw_ij1%has_dijxc_hat==2) then 2822 ii=cplxdij_lmn2_size*ndij 2823 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc_hat,(/ii/)) 2824 indx_dp=indx_dp+ii 2825 end if 2826 if (paw_ij1%has_dijxc_val==2) then 2827 ii=cplxdij_lmn2_size*ndij 2828 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%dijxc_val,(/ii/)) 2829 indx_dp=indx_dp+ii 2830 end if 2831 if (paw_ij1%has_pawu_occ>=1)then 2832 buf_int(indx_int)=size(paw_ij1%noccmmp,1) ;indx_int=indx_int+1 2833 buf_int(indx_int)=size(paw_ij1%noccmmp,2) ;indx_int=indx_int+1 2834 buf_int(indx_int)=size(paw_ij1%noccmmp,3) ;indx_int=indx_int+1 2835 buf_int(indx_int)=size(paw_ij1%noccmmp,4) ;indx_int=indx_int+1 2836 nocc=paw_ij1%ndij 2837 buf_dp(indx_dp:indx_dp+nocc-1)=paw_ij1%nocctot(1:nocc) 2838 indx_dp=indx_dp+nocc 2839 nocc=size(paw_ij1%noccmmp) 2840 buf_dp(indx_dp:indx_dp+nocc-1)=reshape(paw_ij1%noccmmp,(/nocc/)) 2841 indx_dp=indx_dp+nocc 2842 end if 2843 if (paw_ij1%has_exexch_pot>=1) then 2844 sz1=size(paw_ij1%vpawx,1);sz2=size(paw_ij1%vpawx,2) 2845 sz3=size(paw_ij1%vpawx,3) 2846 buf_int(indx_int)=sz1; indx_int=indx_int+1 2847 buf_int(indx_int)=sz2; indx_int=indx_int+1 2848 buf_int(indx_int)=sz3; indx_int=indx_int+1 2849 ii=sz1*sz2*sz3 2850 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij1%vpawx,(/(ii)/)) 2851 indx_dp=indx_dp+ii 2852 end if 2853 end do 2854 indx_int=indx_int-1;indx_dp=indx_dp-1 2855 if ((indx_int/=buf_int_size).or.(indx_dp/=buf_dp_size)) then 2856 write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size 2857 LIBPAW_BUG(msg) 2858 end if 2859 2860 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.
SOURCE
729 subroutine paw_ij_copy(paw_ij_in,paw_ij_cpy, & 730 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 731 732 !Arguments ------------------------------------ 733 !scalars 734 integer,optional,intent(in) :: comm_atom 735 !arrays 736 integer,optional,target,intent(in) :: mpi_atmtab(:) 737 type(paw_ij_type),intent(in) :: paw_ij_in(:) 738 type(paw_ij_type),intent(inout),target :: paw_ij_cpy(:) 739 740 !Local variables------------------------------- 741 !scalars 742 integer :: ij,ij1,my_comm_atom,my_paw_ij,npaw_ij_in,npaw_ij_max,npaw_ij_out,paral_case,sz1,sz2,sz3,sz4 743 logical :: my_atmtab_allocated,paral_atom 744 character(len=500) :: msg 745 !arrays 746 integer,pointer :: my_atmtab(:) 747 type(paw_ij_type),pointer :: paw_ij_out(:) 748 749 ! ************************************************************************* 750 751 !@Paw_ij_type 752 753 !Retrieve sizes 754 npaw_ij_in=size(paw_ij_in);npaw_ij_out=size(paw_ij_cpy) 755 756 !Set up parallelism over atoms 757 paral_atom=(present(comm_atom));if (paral_atom) paral_atom=(xmpi_comm_size(comm_atom)>1) 758 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 759 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 760 my_atmtab_allocated=.false. 761 762 !Determine in which case we are (parallelism, ...) 763 !No parallelism: a single copy operation 764 paral_case=0;npaw_ij_max=npaw_ij_in 765 paw_ij_out => paw_ij_cpy 766 if (paral_atom) then 767 if (npaw_ij_out<npaw_ij_in) then ! Parallelism: the copy operation is a scatter 768 call get_my_natom(my_comm_atom,my_paw_ij,npaw_ij_in) 769 if (my_paw_ij==npaw_ij_out) then 770 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,npaw_ij_in) 771 paral_case=1;npaw_ij_max=npaw_ij_out 772 paw_ij_out => paw_ij_cpy 773 else 774 msg=' npaw_ij_out should be equal to my_paw_ij !' 775 LIBPAW_BUG(msg) 776 end if 777 else ! Parallelism: the copy operation is a gather 778 call get_my_natom(my_comm_atom,my_paw_ij,npaw_ij_out) 779 if (my_paw_ij==npaw_ij_in) then 780 paral_case=2;npaw_ij_max=npaw_ij_in 781 else 782 msg=' npaw_ij_in should be equal to my_paw_ij !' 783 LIBPAW_BUG(msg) 784 end if 785 end if 786 end if 787 788 !First case: a simple copy or a scatter 789 if (npaw_ij_max>0.and.((paral_case==0).or.(paral_case==1))) then 790 call paw_ij_nullify(paw_ij_out) 791 792 ! Loop on paw_ij components 793 do ij1=1,npaw_ij_max 794 ij=ij1; if (paral_case==1) ij=my_atmtab(ij1) 795 796 paw_ij_out(ij1)%qphase=paw_ij_in(ij)%qphase 797 paw_ij_out(ij1)%cplex_dij=paw_ij_in(ij)%cplex_dij 798 paw_ij_out(ij1)%has_dij=paw_ij_in(ij)%has_dij 799 paw_ij_out(ij1)%has_dij0=paw_ij_in(ij)%has_dij0 800 paw_ij_out(ij1)%has_dijexxc=paw_ij_in(ij)%has_dijexxc 801 paw_ij_out(ij1)%has_dijfock=paw_ij_in(ij)%has_dijfock 802 paw_ij_out(ij1)%has_dijfr=paw_ij_in(ij)%has_dijfr 803 paw_ij_out(ij1)%has_dijhartree=paw_ij_in(ij)%has_dijhartree 804 paw_ij_out(ij1)%has_dijhat=paw_ij_in(ij)%has_dijhat 805 paw_ij_out(ij1)%has_dijnd=paw_ij_in(ij)%has_dijnd 806 paw_ij_out(ij1)%has_dijso=paw_ij_in(ij)%has_dijso 807 paw_ij_out(ij1)%has_dijU=paw_ij_in(ij)%has_dijU 808 paw_ij_out(ij1)%has_dijxc=paw_ij_in(ij)%has_dijxc 809 paw_ij_out(ij1)%has_dijxc_hat=paw_ij_in(ij)%has_dijxc_hat 810 paw_ij_out(ij1)%has_dijxc_val=paw_ij_in(ij)%has_dijxc_val 811 paw_ij_out(ij1)%has_exexch_pot=paw_ij_in(ij)%has_exexch_pot 812 paw_ij_out(ij1)%has_pawu_occ=paw_ij_in(ij)%has_pawu_occ 813 paw_ij_out(ij1)%itypat=paw_ij_in(ij)%itypat 814 paw_ij_out(ij1)%lmn_size=paw_ij_in(ij)%lmn_size 815 paw_ij_out(ij1)%lmn2_size=paw_ij_in(ij)%lmn2_size 816 paw_ij_out(ij1)%ndij=paw_ij_in(ij)%ndij 817 paw_ij_out(ij1)%nspden=paw_ij_in(ij)%nspden 818 paw_ij_out(ij1)%nsppol=paw_ij_in(ij)%nsppol 819 if (paw_ij_in(ij)%has_dij>=1) then 820 sz1=size(paw_ij_in(ij)%dij,1);sz2=size(paw_ij_in(ij)%dij,2) 821 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dij,(sz1,sz2)) 822 if (paw_ij_in(ij)%has_dij==2) & 823 & paw_ij_out(ij1)%dij(:,:)=paw_ij_in(ij)%dij(:,:) 824 end if 825 if (paw_ij_in(ij)%has_dij0>=1) then 826 sz1=size(paw_ij_in(ij)%dij0,1) 827 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dij0,(sz1)) 828 if (paw_ij_in(ij)%has_dij0 ==2) & 829 & paw_ij_out(ij1)%dij0(:)=paw_ij_in(ij)%dij0(:) 830 end if 831 if (paw_ij_in(ij)%has_dijexxc>=1) then 832 sz1=size(paw_ij_in(ij)%dijexxc,1);sz2=size(paw_ij_in(ij)%dijexxc,2) 833 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijexxc,(sz1,sz2)) 834 if (paw_ij_in(ij)%has_dijexxc==2) & 835 & paw_ij_out(ij1)%dijexxc(:,:)=paw_ij_in(ij)%dijexxc(:,:) 836 end if 837 if (paw_ij_in(ij)%has_dijfock>=1) then 838 sz1=size(paw_ij_in(ij)%dijfock,1);sz2=size(paw_ij_in(ij)%dijfock,2) 839 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijfock,(sz1,sz2)) 840 if (paw_ij_in(ij)%has_dijfock==2) & 841 & paw_ij_out(ij1)%dijfock(:,:)=paw_ij_in(ij)%dijfock(:,:) 842 end if 843 if (paw_ij_in(ij)%has_dijfr>=1) then 844 sz1=size(paw_ij_in(ij)%dijfr,1);sz2=size(paw_ij_in(ij)%dijfr,2) 845 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijfr,(sz1,sz2)) 846 if (paw_ij_in(ij)%has_dijfr==2) & 847 & paw_ij_out(ij1)%dijfr(:,:)=paw_ij_in(ij)%dijfr(:,:) 848 end if 849 if (paw_ij_in(ij)%has_dijhartree>=1) then 850 sz1=size(paw_ij_in(ij)%dijhartree,1) 851 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijhartree,(sz1)) 852 if (paw_ij_in(ij)%has_dijhartree==2) & 853 & paw_ij_out(ij1)%dijhartree(:)=paw_ij_in(ij)%dijhartree(:) 854 end if 855 if (paw_ij_in(ij)%has_dijhat>=1) then 856 sz1=size(paw_ij_in(ij)%dijhat,1);sz2=size(paw_ij_in(ij)%dijhat,2) 857 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijhat,(sz1,sz2)) 858 if (paw_ij_in(ij)%has_dijhat==2) & 859 & paw_ij_out(ij1)%dijhat(:,:)=paw_ij_in(ij)%dijhat(:,:) 860 end if 861 if (paw_ij_in(ij)%has_dijnd>=1) then 862 sz1=size(paw_ij_in(ij)%dijnd,1);sz2=size(paw_ij_in(ij)%dijnd,2) 863 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijnd,(sz1,sz2)) 864 if (paw_ij_in(ij)%has_dijnd==2) & 865 & paw_ij_out(ij1)%dijnd(:,:)=paw_ij_in(ij)%dijnd(:,:) 866 end if 867 if (paw_ij_in(ij)%has_dijU>=1) then 868 sz1=size(paw_ij_in(ij)%dijU,1);sz2=size(paw_ij_in(ij)%dijU,2) 869 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijU,(sz1,sz2)) 870 if (paw_ij_in(ij)%has_dijU==2) & 871 & paw_ij_out(ij1)%dijU(:,:)=paw_ij_in(ij)%dijU(:,:) 872 end if 873 if (paw_ij_in(ij)%has_dijso>=1) then 874 sz1=size(paw_ij_in(ij)%dijso,1);sz2=size(paw_ij_in(ij)%dijso,2) 875 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijso,(sz1,sz2)) 876 if (paw_ij_in(ij)%has_dijso==2) & 877 & paw_ij_out(ij1)%dijso(:,:)=paw_ij_in(ij)%dijso(:,:) 878 end if 879 if (paw_ij_in(ij)%has_dijxc>=1) then 880 sz1=size(paw_ij_in(ij)%dijxc,1);sz2=size(paw_ij_in(ij)%dijxc,2) 881 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc,(sz1,sz2)) 882 if (paw_ij_in(ij)%has_dijxc==2) & 883 & paw_ij_out(ij1)%dijxc(:,:)=paw_ij_in(ij)%dijxc(:,:) 884 end if 885 if (paw_ij_in(ij)%has_dijxc_hat>=1) then 886 sz1=size(paw_ij_in(ij)%dijxc_hat,1);sz2=size(paw_ij_in(ij)%dijxc_hat,2) 887 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc_hat,(sz1,sz2)) 888 if (paw_ij_in(ij)%has_dijxc_hat==2) & 889 & paw_ij_out(ij1)%dijxc_hat(:,:)=paw_ij_in(ij)%dijxc_hat(:,:) 890 end if 891 if (paw_ij_in(ij)%has_dijxc_val>=1) then 892 sz1=size(paw_ij_in(ij)%dijxc_val,1);sz2=size(paw_ij_in(ij)%dijxc_val,2) 893 LIBPAW_ALLOCATE(paw_ij_out(ij1)%dijxc_val,(sz1,sz2)) 894 if (paw_ij_in(ij)%has_dijxc_val==2) & 895 & paw_ij_out(ij1)%dijxc_val(:,:)=paw_ij_in(ij)%dijxc_val(:,:) 896 end if 897 if (paw_ij_in(ij)%has_pawu_occ>=1) then 898 sz1=size(paw_ij_in(ij)%noccmmp,1);sz2=size(paw_ij_in(ij)%noccmmp,2) 899 sz3=size(paw_ij_in(ij)%noccmmp,3);sz4=size(paw_ij_in(ij)%noccmmp,4) 900 LIBPAW_ALLOCATE(paw_ij_out(ij1)%noccmmp,(sz1,sz2,sz3,sz4)) 901 sz1=size(paw_ij_in(ij)%nocctot,1) 902 LIBPAW_ALLOCATE(paw_ij_out(ij1)%nocctot,(sz1)) 903 if (paw_ij_in(ij)%has_pawu_occ==2) then 904 paw_ij_out(ij1)%noccmmp(:,:,:,:)=paw_ij_in(ij)%noccmmp(:,:,:,:) 905 paw_ij_out(ij1)%nocctot(:)=paw_ij_in(ij)%nocctot(:) 906 end if 907 end if 908 if (paw_ij_in(ij)%has_exexch_pot >= 1) then 909 sz1=size(paw_ij_in(ij)%vpawx,1);sz2=size(paw_ij_in(ij)%vpawx,2) 910 sz3=size(paw_ij_in(ij)%vpawx,3) 911 LIBPAW_ALLOCATE(paw_ij_out(ij1)%vpawx,(sz1,sz2,sz3)) 912 if (paw_ij_in(ij)%has_exexch_pot==2) & 913 & paw_ij_out(ij1)%vpawx(:,:,:)=paw_ij_in(ij)%vpawx(:,:,:) 914 end if 915 916 end do 917 end if 918 919 !Second case: a gather 920 if (paral_case==2) then 921 call paw_ij_gather(paw_ij_in,paw_ij_cpy,-1,my_comm_atom) 922 end if 923 924 !Destroy atom table used for parallelism 925 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 926 927 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
SOURCE
561 subroutine paw_ij_free(Paw_ij) 562 563 !Arguments ------------------------------------ 564 !arrays 565 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 566 567 !Local variables------------------------------- 568 integer :: iat,natom 569 570 ! ************************************************************************* 571 572 !@Paw_ij_type 573 574 natom=SIZE(Paw_ij);if (natom==0) return 575 576 do iat=1,natom 577 if (allocated(Paw_ij(iat)%dij )) then 578 LIBPAW_DEALLOCATE(Paw_ij(iat)%dij) 579 end if 580 if (allocated(Paw_ij(iat)%dij0 )) then 581 LIBPAW_DEALLOCATE(Paw_ij(iat)%dij0) 582 end if 583 if (allocated(Paw_ij(iat)%dijexxc )) then 584 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijexxc) 585 end if 586 if (allocated(Paw_ij(iat)%dijfock )) then 587 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijfock) 588 end if 589 if (allocated(Paw_ij(iat)%dijfr )) then 590 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijfr) 591 end if 592 if (allocated(Paw_ij(iat)%dijhartree)) then 593 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijhartree) 594 end if 595 if (allocated(Paw_ij(iat)%dijhat )) then 596 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijhat) 597 end if 598 if (allocated(Paw_ij(iat)%dijnd )) then 599 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijnd) 600 end if 601 if (allocated(Paw_ij(iat)%dijU )) then 602 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijU) 603 end if 604 if (allocated(Paw_ij(iat)%dijso )) then 605 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijso) 606 end if 607 if (allocated(Paw_ij(iat)%dijxc )) then 608 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc) 609 end if 610 if (allocated(Paw_ij(iat)%dijxc_hat )) then 611 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc_hat) 612 end if 613 if (allocated(Paw_ij(iat)%dijxc_val )) then 614 LIBPAW_DEALLOCATE(Paw_ij(iat)%dijxc_val) 615 end if 616 if (allocated(Paw_ij(iat)%noccmmp )) then 617 LIBPAW_DEALLOCATE(Paw_ij(iat)%noccmmp) 618 end if 619 if (allocated(Paw_ij(iat)%nocctot )) then 620 LIBPAW_DEALLOCATE(Paw_ij(iat)%nocctot) 621 end if 622 if (allocated(Paw_ij(iat)%vpawx )) then 623 LIBPAW_DEALLOCATE(Paw_ij(iat)%vpawx) 624 end if 625 626 ! === Reset all has_* flags === 627 Paw_ij(iat)%has_dij =0 628 Paw_ij(iat)%has_dij0 =0 629 Paw_ij(iat)%has_dijexxc =0 630 Paw_ij(iat)%has_dijfock =0 631 Paw_ij(iat)%has_dijfr =0 632 Paw_ij(iat)%has_dijhartree=0 633 Paw_ij(iat)%has_dijhat =0 634 Paw_ij(iat)%has_dijnd =0 635 Paw_ij(iat)%has_dijso =0 636 Paw_ij(iat)%has_dijU =0 637 Paw_ij(iat)%has_dijxc =0 638 Paw_ij(iat)%has_dijxc_hat =0 639 Paw_ij(iat)%has_dijxc_val =0 640 Paw_ij(iat)%has_exexch_pot=0 641 Paw_ij(iat)%has_pawu_occ =0 642 end do 643 644 !call paw_ij_nullify(Paw_ij) 645 646 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
SOURCE
1324 subroutine paw_ij_gather(paw_ij_in,paw_ij_gathered,master,comm_atom) 1325 1326 !Arguments ------------------------------------ 1327 !scalars 1328 integer,intent(in) :: master,comm_atom 1329 !arrays 1330 type(paw_ij_type),intent(in) :: paw_ij_in(:) 1331 type(paw_ij_type),intent(inout) :: paw_ij_gathered(:) 1332 1333 !Local variables------------------------------- 1334 !scalars 1335 integer :: buf_dp_size,buf_dp_size_all,buf_int_size,buf_int_size_all 1336 integer :: cplxdij_lmn2_size,cplxdijq_lmn2_size,cplxq_lmn2_size 1337 integer :: iat,ii,ierr,ij,indx_dp,indx_int,lmn2_size 1338 integer :: me_atom,ndij,nocc,nocc1,nocc2,nocc3,nocc4,npaw_ij_in,npaw_ij_in_sum 1339 integer :: npaw_ij_out,nproc_atom,nspden,sz1,sz2,sz3,sz4 1340 logical :: my_atmtab_allocated,paral_atom 1341 character(len=500) :: msg 1342 !arrays 1343 integer :: bufsz(2) 1344 integer,allocatable :: buf_int(:),buf_int_all(:) 1345 integer,allocatable :: count_dp(:),count_int(:),count_tot(:),displ_dp(:),displ_int(:) 1346 integer,pointer :: my_atmtab(:) 1347 real(dp),allocatable :: buf_dp(:),buf_dp_all(:) 1348 1349 ! ************************************************************************* 1350 1351 !@Paw_ij_type 1352 1353 npaw_ij_in=size(paw_ij_in);npaw_ij_out=size(paw_ij_gathered) 1354 1355 nproc_atom=xmpi_comm_size(comm_atom) 1356 me_atom=xmpi_comm_rank(comm_atom) 1357 1358 !Special case 1 process 1359 if (nproc_atom==1) then 1360 if (master==-1.or.me_atom==master) then 1361 call paw_ij_free(paw_ij_gathered) 1362 call paw_ij_nullify(paw_ij_gathered) 1363 do iat=1,npaw_ij_in 1364 paw_ij_gathered(iat)%cplex_dij =paw_ij_in(iat)%cplex_dij 1365 paw_ij_gathered(iat)%qphase =paw_ij_in(iat)%qphase 1366 Paw_ij_gathered(iat)%has_dij =paw_ij_in(iat)%has_dij 1367 Paw_ij_gathered(iat)%has_dij0 =paw_ij_in(iat)%has_dij0 1368 Paw_ij_gathered(iat)%has_dijexxc =paw_ij_in(iat)%has_dijexxc 1369 Paw_ij_gathered(iat)%has_dijfock =paw_ij_in(iat)%has_dijfock 1370 Paw_ij_gathered(iat)%has_dijfr =paw_ij_in(iat)%has_dijfr 1371 Paw_ij_gathered(iat)%has_dijhartree=paw_ij_in(iat)%has_dijhartree 1372 Paw_ij_gathered(iat)%has_dijhat =paw_ij_in(iat)%has_dijhat 1373 Paw_ij_gathered(iat)%has_dijnd =paw_ij_in(iat)%has_dijnd 1374 Paw_ij_gathered(iat)%has_dijso =paw_ij_in(iat)%has_dijso 1375 Paw_ij_gathered(iat)%has_dijU =paw_ij_in(iat)%has_dijU 1376 Paw_ij_gathered(iat)%has_dijxc =paw_ij_in(iat)%has_dijxc 1377 Paw_ij_gathered(iat)%has_dijxc_hat =paw_ij_in(iat)%has_dijxc_hat 1378 Paw_ij_gathered(iat)%has_dijxc_val =paw_ij_in(iat)%has_dijxc_val 1379 Paw_ij_gathered(iat)%has_exexch_pot=paw_ij_in(iat)%has_exexch_pot 1380 Paw_ij_gathered(iat)%has_pawu_occ =paw_ij_in(iat)%has_pawu_occ 1381 paw_ij_gathered(iat)%itypat =paw_ij_in(iat)%itypat 1382 paw_ij_gathered(iat)%lmn_size =paw_ij_in(iat)%lmn_size 1383 paw_ij_gathered(iat)%lmn2_size =paw_ij_in(iat)%lmn2_size 1384 paw_ij_gathered(iat)%ndij =paw_ij_in(iat)%ndij 1385 paw_ij_gathered(iat)%nspden =paw_ij_in(iat)%nspden 1386 paw_ij_gathered(iat)%nsppol =paw_ij_in(iat)%nsppol 1387 lmn2_size=paw_ij_gathered(iat)%lmn2_size 1388 cplxdij_lmn2_size=paw_ij_gathered(iat)%cplex_dij*lmn2_size 1389 cplxq_lmn2_size=paw_ij_gathered(iat)%qphase*lmn2_size 1390 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij_gathered(iat)%qphase 1391 ndij=paw_ij_gathered(iat)%ndij 1392 nspden=paw_ij_gathered(iat)%nspden 1393 if (paw_ij_gathered(iat)%has_dij>=1) then 1394 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij,(cplxdijq_lmn2_size,ndij)) 1395 if (paw_ij_in(iat)%has_dij==2) then 1396 paw_ij_gathered(iat)%dij=paw_ij_in(iat)%dij 1397 end if 1398 end if 1399 if (paw_ij_gathered(iat)%has_dij0 >=1) then 1400 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij0,(lmn2_size)) 1401 if (paw_ij_in(iat)%has_dij0==2) then 1402 paw_ij_gathered(iat)%dij0=paw_ij_in(iat)%dij0 1403 end if 1404 end if 1405 if (paw_ij_gathered(iat)%has_dijexxc >=1) then 1406 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijexxc,(cplxdij_lmn2_size,ndij)) 1407 if (paw_ij_in(iat)%has_dijexxc==2) then 1408 paw_ij_gathered(iat)%dijexxc(:,:)=paw_ij_in(iat)%dijexxc(:,:) 1409 end if 1410 end if 1411 if (paw_ij_gathered(iat)%has_dijfock >=1) then 1412 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfock,(cplxdij_lmn2_size,ndij)) 1413 if (paw_ij_in(iat)%has_dijfock==2) then 1414 paw_ij_gathered(iat)%dijfock(:,:)=paw_ij_in(iat)%dijfock(:,:) 1415 end if 1416 end if 1417 if (paw_ij_gathered(iat)%has_dijfr >=1) then 1418 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfr,(cplxdijq_lmn2_size,ndij)) 1419 if (paw_ij_in(iat)%has_dijfr==2) then 1420 paw_ij_gathered(iat)%dijfr(:,:)=paw_ij_in(iat)%dijfr(:,:) 1421 end if 1422 end if 1423 if (paw_ij_gathered(iat)%has_dijhartree >=1) then 1424 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhartree,(cplxq_lmn2_size)) 1425 if (paw_ij_in(iat)%has_dijhartree==2) then 1426 paw_ij_gathered(iat)%dijhartree(:)=paw_ij_in(iat)%dijhartree(:) 1427 end if 1428 end if 1429 if (paw_ij_gathered(iat)%has_dijhat >=1) then 1430 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhat,(cplxdijq_lmn2_size,ndij)) 1431 if (paw_ij_in(iat)%has_dijhat==2) then 1432 paw_ij_gathered(iat)%dijhat(:,:)=paw_ij_in(iat)%dijhat(:,:) 1433 end if 1434 end if 1435 if (paw_ij_gathered(iat)%has_dijnd >=1) then 1436 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijnd,(cplxdij_lmn2_size,ndij)) 1437 if (paw_ij_in(iat)%has_dijnd==2) then 1438 paw_ij_gathered(iat)%dijnd(:,:)=paw_ij_in(iat)%dijnd(:,:) 1439 end if 1440 end if 1441 if (paw_ij_gathered(iat)%has_dijU >=1) then 1442 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijU,(cplxdijq_lmn2_size,ndij)) 1443 if (paw_ij_in(iat)%has_dijU==2) then 1444 paw_ij_gathered(iat)%dijU(:,:)=paw_ij_in(iat)%dijU(:,:) 1445 end if 1446 end if 1447 if (paw_ij_gathered(iat)%has_dijso >=1) then 1448 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijso,(cplxdijq_lmn2_size,ndij)) 1449 if (paw_ij_in(iat)%has_dijso==2) then 1450 paw_ij_gathered(iat)%dijso(:,:)=paw_ij_in(iat)%dijso(:,:) 1451 end if 1452 end if 1453 if (paw_ij_gathered(iat)%has_dijxc >=1) then 1454 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc,(cplxdijq_lmn2_size,ndij)) 1455 if (paw_ij_in(iat)%has_dijxc==2) then 1456 paw_ij_gathered(iat)%dijxc(:,:)=paw_ij_in(iat)%dijxc(:,:) 1457 end if 1458 end if 1459 if (paw_ij_gathered(iat)%has_dijxc_hat >=1) then 1460 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_hat,(cplxdij_lmn2_size,ndij)) 1461 if (paw_ij_in(iat)%has_dijxc_hat==2) then 1462 paw_ij_gathered(iat)%dijxc_hat(:,:)=paw_ij_in(iat)%dijxc_hat(:,:) 1463 end if 1464 end if 1465 if (paw_ij_gathered(iat)%has_dijxc_val >=1) then 1466 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_val,(cplxdij_lmn2_size,ndij)) 1467 if (paw_ij_in(iat)%has_dijxc_val==2) then 1468 paw_ij_gathered(iat)%dijxc_val(:,:)=paw_ij_in(iat)%dijxc_val(:,:) 1469 end if 1470 end if 1471 if (paw_ij_gathered(iat)%has_pawu_occ >=1) then 1472 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%nocctot,(ndij)) 1473 if (paw_ij_gathered(iat)%has_pawu_occ==2) then 1474 paw_ij_gathered(iat)%nocctot(:)=paw_ij_in(iat)%nocctot(:) 1475 if (allocated(paw_ij_in(iat)%noccmmp)) then 1476 sz1=size(paw_ij_in(iat)%noccmmp,1);sz2=size(paw_ij_in(iat)%noccmmp,2) 1477 sz3=size(paw_ij_in(iat)%noccmmp,3);sz4=size(paw_ij_in(iat)%noccmmp,4) 1478 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%noccmmp,(sz1,sz2,sz3,sz4)) 1479 paw_ij_gathered(iat)%noccmmp(:,:,:,:)=paw_ij_in(iat)%noccmmp(:,:,:,:) 1480 end if 1481 end if 1482 end if 1483 if (paw_ij_in(iat)%has_exexch_pot >=1) then 1484 sz1=size(paw_ij_in(iat)%vpawx,1);sz2=size(paw_ij_in(iat)%vpawx,2) 1485 sz3=size(paw_ij_in(iat)%vpawx,3) 1486 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%vpawx,(sz1,sz2,sz3)) 1487 if (paw_ij_in(iat)%has_exexch_pot==2) then 1488 paw_ij_gathered(iat)%vpawx(:,:,:)=paw_ij_in(iat)%vpawx(:,:,:) 1489 end if 1490 end if 1491 end do ! iat 1492 end if 1493 return 1494 end if !nproc_atom =1 1495 1496 !Test on sizes 1497 npaw_ij_in_sum=npaw_ij_in 1498 1499 call xmpi_sum(npaw_ij_in_sum,comm_atom,ierr) 1500 if (master==-1) then 1501 if (npaw_ij_out/=npaw_ij_in_sum) then 1502 msg='Wrong sizes sum[npaw_ij_ij]/=npaw_ij_out !' 1503 LIBPAW_BUG(msg) 1504 end if 1505 else 1506 if (me_atom==master.and.npaw_ij_out/=npaw_ij_in_sum) then 1507 msg='(2) paw_ij_gathered wrongly allocated !' 1508 LIBPAW_BUG(msg) 1509 end if 1510 end if 1511 1512 !Retrieve table of atoms 1513 paral_atom=.true.;nullify(my_atmtab) 1514 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,npaw_ij_in_sum) 1515 1516 !Compute sizes of buffers 1517 buf_int_size=0;buf_dp_size=0 1518 do ij=1,npaw_ij_in 1519 lmn2_size=paw_ij_in(ij)%lmn2_size 1520 cplxdij_lmn2_size=paw_ij_in(ij)%cplex_dij*lmn2_size 1521 cplxq_lmn2_size=paw_ij_in(ij)%qphase*lmn2_size 1522 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij_in(ij)%qphase 1523 ndij=paw_ij_in(ij)%ndij 1524 buf_int_size=buf_int_size+24 1525 if (paw_ij_in(ij)%has_dij==2) then 1526 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1527 end if 1528 if (paw_ij_in(ij)%has_dij0==2) then 1529 buf_dp_size=buf_dp_size +lmn2_size 1530 end if 1531 if (paw_ij_in(ij)%has_dijexxc==2) then 1532 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1533 end if 1534 if (paw_ij_in(ij)%has_dijfock==2) then 1535 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1536 end if 1537 if (paw_ij_in(ij)%has_dijfr==2) then 1538 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1539 end if 1540 if (paw_ij_in(ij)%has_dijhartree==2) then 1541 buf_dp_size=buf_dp_size +cplxq_lmn2_size 1542 end if 1543 if (paw_ij_in(ij)%has_dijhat==2) then 1544 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1545 end if 1546 if (paw_ij_in(ij)%has_dijnd==2) then 1547 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1548 end if 1549 if (paw_ij_in(ij)%has_dijso==2) then 1550 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1551 end if 1552 if (paw_ij_in(ij)%has_dijU==2) then 1553 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1554 end if 1555 if (paw_ij_in(ij)%has_dijxc==2) then 1556 buf_dp_size=buf_dp_size +cplxdijq_lmn2_size*ndij 1557 end if 1558 if (paw_ij_in(ij)%has_dijxc_hat==2) then 1559 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1560 end if 1561 if (paw_ij_in(ij)%has_dijxc_val==2) then 1562 buf_dp_size=buf_dp_size +cplxdij_lmn2_size*ndij 1563 end if 1564 if (paw_ij_in(ij)%has_pawu_occ>=1) then 1565 buf_int_size=buf_int_size+4 1566 buf_dp_size=buf_dp_size & 1567 & +size(paw_ij_in(ij)%nocctot) & 1568 & +size(paw_ij_in(ij)%noccmmp) 1569 end if 1570 if (paw_ij_in(ij)%has_exexch_pot>=1) then 1571 buf_int_size=buf_int_size+3 1572 buf_dp_size=buf_dp_size +size(paw_ij_in(ij)%vpawx) 1573 end if 1574 end do 1575 1576 !Fill input buffers 1577 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 1578 LIBPAW_ALLOCATE(buf_dp ,(buf_dp_size)) 1579 indx_int=1;indx_dp=1 1580 do ij=1,npaw_ij_in 1581 nspden=paw_ij_in(ij)%nspden 1582 buf_int(indx_int)=my_atmtab(ij) ;indx_int=indx_int+1 1583 buf_int(indx_int)=paw_ij_in(ij)%cplex_dij ;indx_int=indx_int+1 1584 buf_int(indx_int)=paw_ij_in(ij)%qphase ;indx_int=indx_int+1 1585 buf_int(indx_int)=paw_ij_in(ij)%itypat ;indx_int=indx_int+1 1586 buf_int(indx_int)=paw_ij_in(ij)%nspden ;indx_int=indx_int+1 1587 buf_int(indx_int)=paw_ij_in(ij)%nsppol ;indx_int=indx_int+1 1588 buf_int(indx_int)=paw_ij_in(ij)%lmn_size ;indx_int=indx_int+1 1589 buf_int(indx_int)=paw_ij_in(ij)%lmn2_size ;indx_int=indx_int+1 1590 buf_int(indx_int)=paw_ij_in(ij)%ndij ;indx_int=indx_int+1 1591 buf_int(indx_int)=paw_ij_in(ij)%has_dij ;indx_int=indx_int+1 1592 buf_int(indx_int)=paw_ij_in(ij)%has_dij0 ;indx_int=indx_int+1 1593 buf_int(indx_int)=paw_ij_in(ij)%has_dijexxc ;indx_int=indx_int+1 1594 buf_int(indx_int)=paw_ij_in(ij)%has_dijfock ;indx_int=indx_int+1 1595 buf_int(indx_int)=paw_ij_in(ij)%has_dijfr ;indx_int=indx_int+1 1596 buf_int(indx_int)=paw_ij_in(ij)%has_dijhartree ;indx_int=indx_int+1 1597 buf_int(indx_int)=paw_ij_in(ij)%has_dijhat ;indx_int=indx_int+1 1598 buf_int(indx_int)=paw_ij_in(ij)%has_dijnd ;indx_int=indx_int+1 1599 buf_int(indx_int)=paw_ij_in(ij)%has_dijso ;indx_int=indx_int+1 1600 buf_int(indx_int)=paw_ij_in(ij)%has_dijU ;indx_int=indx_int+1 1601 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc ;indx_int=indx_int+1 1602 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc_hat ;indx_int=indx_int+1 1603 buf_int(indx_int)=paw_ij_in(ij)%has_dijxc_val ;indx_int=indx_int+1 1604 buf_int(indx_int)=paw_ij_in(ij)%has_exexch_pot ;indx_int=indx_int+1 1605 buf_int(indx_int)=paw_ij_in(ij)%has_pawu_occ ;indx_int=indx_int+1 1606 lmn2_size=paw_ij_in(ij)%lmn2_size 1607 cplxdij_lmn2_size=paw_ij_in(ij)%cplex_dij*lmn2_size 1608 cplxq_lmn2_size=paw_ij_in(ij)%qphase*lmn2_size 1609 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij_in(ij)%qphase 1610 ndij=paw_ij_in(ij)%ndij 1611 if (paw_ij_in(ij)%has_dij==2) then 1612 ii=cplxdijq_lmn2_size*ndij 1613 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dij,(/ii/)) 1614 indx_dp=indx_dp+ii 1615 end if 1616 if (paw_ij_in(ij)%has_dij0==2) then 1617 ii=lmn2_size 1618 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij_in(ij)%dij0(1:ii) 1619 indx_dp=indx_dp+ii 1620 end if 1621 if (paw_ij_in(ij)%has_dijexxc==2) then 1622 ii=cplxdij_lmn2_size*ndij 1623 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijexxc,(/ii/)) 1624 indx_dp=indx_dp+ii 1625 end if 1626 if (paw_ij_in(ij)%has_dijfock==2) then 1627 ii=cplxdij_lmn2_size*ndij 1628 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijfock,(/ii/)) 1629 indx_dp=indx_dp+ii 1630 end if 1631 if (paw_ij_in(ij)%has_dijfr==2) then 1632 ii=cplxdijq_lmn2_size*ndij 1633 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijfr,(/ii/)) 1634 indx_dp=indx_dp+ii 1635 end if 1636 if (paw_ij_in(ij)%has_dijhartree==2) then 1637 ii=cplxq_lmn2_size 1638 buf_dp(indx_dp:indx_dp+ii-1)=paw_ij_in(ij)%dijhartree(1:ii) 1639 indx_dp=indx_dp+ii 1640 end if 1641 if (paw_ij_in(ij)%has_dijhat==2) then 1642 ii=cplxdijq_lmn2_size*ndij 1643 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijhat,(/ii/)) 1644 indx_dp=indx_dp+ii 1645 end if 1646 if (paw_ij_in(ij)%has_dijnd==2) then 1647 ii=cplxdij_lmn2_size*ndij 1648 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijnd,(/ii/)) 1649 indx_dp=indx_dp+ii 1650 end if 1651 if (paw_ij_in(ij)%has_dijso==2) then 1652 ii=cplxdijq_lmn2_size*ndij 1653 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijso,(/ii/)) 1654 indx_dp=indx_dp+ii 1655 end if 1656 if (paw_ij_in(ij)%has_dijU==2) then 1657 ii=cplxdijq_lmn2_size*ndij 1658 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijU,(/ii/)) 1659 indx_dp=indx_dp+ii 1660 end if 1661 if (paw_ij_in(ij)%has_dijxc==2) then 1662 ii=cplxdijq_lmn2_size*ndij 1663 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc,(/ii/)) 1664 indx_dp=indx_dp+ii 1665 end if 1666 if (paw_ij_in(ij)%has_dijxc_hat==2) then 1667 ii=cplxdij_lmn2_size*ndij 1668 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc_hat,(/ii/)) 1669 indx_dp=indx_dp+ii 1670 end if 1671 if (paw_ij_in(ij)%has_dijxc_val==2) then 1672 ii=cplxdij_lmn2_size*ndij 1673 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%dijxc_val,(/ii/)) 1674 indx_dp=indx_dp+ii 1675 end if 1676 if (paw_ij_in(ij)%has_pawu_occ>=1)then 1677 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,1) ;indx_int=indx_int+1 1678 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,2) ;indx_int=indx_int+1 1679 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,3) ;indx_int=indx_int+1 1680 buf_int(indx_int)=size(paw_ij_in(ij)%noccmmp,4) ;indx_int=indx_int+1 1681 nocc=paw_ij_in(ij)%ndij 1682 buf_dp(indx_dp:indx_dp+nocc-1)=paw_ij_in(ij)%nocctot(1:nocc) 1683 indx_dp=indx_dp+nocc 1684 nocc=size(paw_ij_in(ij)%noccmmp) 1685 buf_dp(indx_dp:indx_dp+nocc-1)=reshape(paw_ij_in(ij)%noccmmp,(/nocc/)) 1686 indx_dp=indx_dp+nocc 1687 end if 1688 if (paw_ij_in(ij)%has_exexch_pot>=1) then 1689 sz1=size(paw_ij_in(ij)%vpawx,1);sz2=size(paw_ij_in(ij)%vpawx,2) 1690 sz3=size(paw_ij_in(ij)%vpawx,3) 1691 buf_int(indx_int)=sz1; indx_int=indx_int+1 1692 buf_int(indx_int)=sz2; indx_int=indx_int+1 1693 buf_int(indx_int)=sz3; indx_int=indx_int+1 1694 ii=sz1*sz2*sz3 1695 buf_dp(indx_dp:indx_dp+ii-1)=reshape(paw_ij_in(ij)%vpawx,(/(ii)/)) 1696 indx_dp=indx_dp+ii 1697 end if 1698 end do 1699 1700 !Check 1701 indx_int=indx_int-1;indx_dp=indx_dp-1 1702 if ((indx_int/=buf_int_size).or.(indx_dp/=buf_dp_size)) then 1703 write(msg,*) 'Wrong buffer sizes: buf_int_size=',buf_int_size, & 1704 ' indx_int_size=',indx_int,' buf_dp_size=',buf_dp_size , ' indx_dp=', indx_dp 1705 LIBPAW_BUG(msg) 1706 end if 1707 1708 !Communicate (1 gather for integers, 1 gather for reals) 1709 LIBPAW_ALLOCATE(count_int,(nproc_atom)) 1710 LIBPAW_ALLOCATE(displ_int,(nproc_atom)) 1711 LIBPAW_ALLOCATE(count_dp ,(nproc_atom)) 1712 LIBPAW_ALLOCATE(displ_dp ,(nproc_atom)) 1713 LIBPAW_ALLOCATE(count_tot,(2*nproc_atom)) 1714 bufsz(1)=buf_int_size; bufsz(2)=buf_dp_size 1715 call xmpi_allgather(bufsz,2,count_tot,comm_atom,ierr) 1716 do ij=1,nproc_atom 1717 count_int(ij)=count_tot(2*ij-1) 1718 count_dp (ij)=count_tot(2*ij) 1719 end do 1720 displ_int(1)=0;displ_dp(1)=0 1721 do ij=2,nproc_atom 1722 displ_int(ij)=displ_int(ij-1)+count_int(ij-1) 1723 displ_dp (ij)=displ_dp (ij-1)+count_dp (ij-1) 1724 end do 1725 buf_int_size_all=sum(count_int) 1726 buf_dp_size_all =sum(count_dp) 1727 LIBPAW_DEALLOCATE(count_tot) 1728 if (master==-1.or.me_atom==master) then 1729 LIBPAW_ALLOCATE(buf_int_all,(buf_int_size_all)) 1730 LIBPAW_ALLOCATE(buf_dp_all ,(buf_dp_size_all)) 1731 else 1732 LIBPAW_ALLOCATE(buf_int_all,(0)) 1733 LIBPAW_ALLOCATE(buf_dp_all ,(0)) 1734 end if 1735 if (master==-1) then 1736 call xmpi_allgatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,comm_atom,ierr) 1737 call xmpi_allgatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,comm_atom,ierr) 1738 else 1739 call xmpi_gatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,master,comm_atom,ierr) 1740 call xmpi_gatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,master,comm_atom,ierr) 1741 end if 1742 LIBPAW_DEALLOCATE(count_int) 1743 LIBPAW_DEALLOCATE(displ_int) 1744 LIBPAW_DEALLOCATE(count_dp) 1745 LIBPAW_DEALLOCATE(displ_dp) 1746 1747 !Retrieve data from output buffer 1748 if (master==-1.or.me_atom==master) then 1749 indx_int=1;indx_dp=1 1750 call paw_ij_free(paw_ij_gathered) 1751 call paw_ij_nullify(paw_ij_gathered) 1752 do ij=1,npaw_ij_out 1753 iat=buf_int_all(indx_int) ;indx_int=indx_int+1 1754 paw_ij_gathered(iat)%cplex_dij=buf_int_all(indx_int) ;indx_int=indx_int+1 1755 paw_ij_gathered(iat)%qphase=buf_int_all(indx_int) ;indx_int=indx_int+1 1756 paw_ij_gathered(iat)%itypat=buf_int_all(indx_int) ;indx_int=indx_int+1 1757 paw_ij_gathered(iat)%nspden=buf_int_all(indx_int) ;indx_int=indx_int+1 1758 paw_ij_gathered(iat)%nsppol=buf_int_all(indx_int) ;indx_int=indx_int+1 1759 paw_ij_gathered(iat)%lmn_size=buf_int_all(indx_int) ;indx_int=indx_int+1 1760 paw_ij_gathered(iat)%lmn2_size=buf_int_all(indx_int) ;indx_int=indx_int+1 1761 paw_ij_gathered(iat)%ndij=buf_int_all(indx_int) ;indx_int=indx_int+1 1762 paw_ij_gathered(iat)%has_dij=buf_int_all(indx_int) ;indx_int=indx_int+1 1763 paw_ij_gathered(iat)%has_dij0=buf_int_all(indx_int) ;indx_int=indx_int+1 1764 paw_ij_gathered(iat)%has_dijexxc=buf_int_all(indx_int) ;indx_int=indx_int+1 1765 paw_ij_gathered(iat)%has_dijfock=buf_int_all(indx_int) ;indx_int=indx_int+1 1766 paw_ij_gathered(iat)%has_dijfr=buf_int_all(indx_int) ;indx_int=indx_int+1 1767 paw_ij_gathered(iat)%has_dijhartree=buf_int_all(indx_int) ;indx_int=indx_int+1 1768 paw_ij_gathered(iat)%has_dijhat=buf_int_all(indx_int) ;indx_int=indx_int+1 1769 paw_ij_gathered(iat)%has_dijnd=buf_int_all(indx_int) ;indx_int=indx_int+1 1770 paw_ij_gathered(iat)%has_dijso=buf_int_all(indx_int) ;indx_int=indx_int+1 1771 paw_ij_gathered(iat)%has_dijU=buf_int_all(indx_int) ;indx_int=indx_int+1 1772 paw_ij_gathered(iat)%has_dijxc=buf_int_all(indx_int) ;indx_int=indx_int+1 1773 paw_ij_gathered(iat)%has_dijxc_hat=buf_int_all(indx_int) ;indx_int=indx_int+1 1774 paw_ij_gathered(iat)%has_dijxc_val=buf_int_all(indx_int) ;indx_int=indx_int+1 1775 paw_ij_gathered(iat)%has_exexch_pot=buf_int_all(indx_int) ;indx_int=indx_int+1 1776 paw_ij_gathered(iat)%has_pawu_occ=buf_int_all(indx_int) ;indx_int=indx_int+1 1777 if (paw_ij_gathered(iat)%has_pawu_occ>=1) then 1778 nocc1=buf_int_all(indx_int) ;indx_int=indx_int+1 1779 nocc2=buf_int_all(indx_int) ;indx_int=indx_int+1 1780 nocc3=buf_int_all(indx_int) ;indx_int=indx_int+1 1781 nocc4=buf_int_all(indx_int) ;indx_int=indx_int+1 1782 else 1783 nocc1=0;nocc2=0;nocc3=0;nocc4=0 1784 end if 1785 lmn2_size=paw_ij_gathered(iat)%lmn2_size 1786 cplxdij_lmn2_size=paw_ij_gathered(iat)%cplex_dij*lmn2_size 1787 cplxq_lmn2_size=paw_ij_gathered(iat)%qphase*lmn2_size 1788 cplxdijq_lmn2_size=cplxdij_lmn2_size*paw_ij_gathered(iat)%qphase 1789 ndij=paw_ij_gathered(iat)%ndij 1790 1791 if (paw_ij_gathered(iat)%has_dij>=1) then 1792 ii=cplxdijq_lmn2_size 1793 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij,(ii,ndij)) 1794 if (paw_ij_gathered(iat)%has_dij==2) then 1795 paw_ij_gathered(iat)%dij(:,:)= & 1796 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1797 indx_dp=indx_dp+ii*ndij 1798 end if 1799 end if 1800 if (paw_ij_gathered(iat)%has_dij0 >=1) then 1801 ii=lmn2_size 1802 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dij0,(ii)) 1803 if (paw_ij_gathered(iat)%has_dij0==2) then 1804 paw_ij_gathered(iat)%dij0(:)=buf_dp_all(indx_dp:indx_dp+ii-1) 1805 indx_dp=indx_dp+ii 1806 end if 1807 end if 1808 if (paw_ij_gathered(iat)%has_dijexxc >=1) then 1809 ii=cplxdij_lmn2_size 1810 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijexxc,(ii,ndij)) 1811 if (paw_ij_gathered(iat)%has_dijexxc==2) then 1812 paw_ij_gathered(iat)%dijexxc(:,:)= & 1813 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1814 indx_dp=indx_dp+ii*ndij 1815 end if 1816 end if 1817 if (paw_ij_gathered(iat)%has_dijfock >=1) then 1818 ii=cplxdij_lmn2_size 1819 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfock,(ii,ndij)) 1820 if (paw_ij_gathered(iat)%has_dijfock==2) then 1821 paw_ij_gathered(iat)%dijfock(:,:)= & 1822 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1823 indx_dp=indx_dp+ii*ndij 1824 end if 1825 end if 1826 if (paw_ij_gathered(iat)%has_dijfr >=1) then 1827 ii=cplxdijq_lmn2_size 1828 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijfr,(ii,ndij)) 1829 if (paw_ij_gathered(iat)%has_dijfr==2) then 1830 paw_ij_gathered(iat)%dijfr(:,:)= & 1831 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1832 indx_dp=indx_dp+ii*ndij 1833 end if 1834 end if 1835 if (paw_ij_gathered(iat)%has_dijhartree >=1) then 1836 ii=cplxq_lmn2_size 1837 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhartree,(ii)) 1838 if (paw_ij_gathered(iat)%has_dijhartree==2) then 1839 paw_ij_gathered(iat)%dijhartree(:)=buf_dp_all(indx_dp:indx_dp+ii-1) 1840 indx_dp=indx_dp+ii 1841 end if 1842 end if 1843 if (paw_ij_gathered(iat)%has_dijhat >=1) then 1844 ii=cplxdijq_lmn2_size 1845 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijhat,(ii,ndij)) 1846 if (paw_ij_gathered(iat)%has_dijhat==2) then 1847 paw_ij_gathered(iat)%dijhat(:,:)= & 1848 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1849 indx_dp=indx_dp+ii*ndij 1850 end if 1851 endif 1852 if (paw_ij_gathered(iat)%has_dijnd >=1) then 1853 ii=cplxdij_lmn2_size 1854 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijnd,(ii,ndij)) 1855 if (paw_ij_gathered(iat)%has_dijnd==2) then 1856 paw_ij_gathered(iat)%dijnd(:,:)= & 1857 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1858 indx_dp=indx_dp+ii*ndij 1859 end if 1860 end if 1861 if (paw_ij_gathered(iat)%has_dijso >=1) then 1862 ii=cplxdijq_lmn2_size 1863 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijso,(ii,ndij)) 1864 if (paw_ij_gathered(iat)%has_dijso==2) then 1865 paw_ij_gathered(iat)%dijso(:,:)= & 1866 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1867 indx_dp=indx_dp+ii*ndij 1868 end if 1869 end if 1870 if (paw_ij_gathered(iat)%has_dijU >=1) then 1871 ii=cplxdijq_lmn2_size 1872 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijU,(ii,ndij)) 1873 if (paw_ij_gathered(iat)%has_dijU==2) then 1874 paw_ij_gathered(iat)%dijU(:,:)= & 1875 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1876 indx_dp=indx_dp+ii*ndij 1877 end if 1878 end if 1879 if (paw_ij_gathered(iat)%has_dijxc >=1) then 1880 ii=cplxdijq_lmn2_size 1881 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc,(ii,ndij)) 1882 if (paw_ij_gathered(iat)%has_dijxc==2) then 1883 paw_ij_gathered(iat)%dijxc(:,:)= & 1884 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1885 indx_dp=indx_dp+ii*ndij 1886 end if 1887 end if 1888 if (paw_ij_gathered(iat)%has_dijxc_hat >=1) then 1889 ii=cplxdij_lmn2_size 1890 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_hat,(ii,ndij)) 1891 if (paw_ij_gathered(iat)%has_dijxc_hat==2) then 1892 paw_ij_gathered(iat)%dijxc_hat(:,:)= & 1893 & reshape(buf_dp_all(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 1894 indx_dp=indx_dp+ii*ndij 1895 end if 1896 end if 1897 if (paw_ij_gathered(iat)%has_dijxc_val >=1) then 1898 ii=cplxdij_lmn2_size 1899 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%dijxc_val,(ii,ndij)) 1900 if (paw_ij_gathered(iat)%has_dijxc_val==2) then 1901 paw_ij_gathered(iat)%dijxc_val(:,:)= & 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_pawu_occ >=1) then 1907 nocc=paw_ij_gathered(iat)%ndij 1908 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%nocctot,(nocc)) 1909 paw_ij_gathered(iat)%nocctot(1:nocc)=buf_dp_all(indx_dp:indx_dp+nocc-1) 1910 indx_dp=indx_dp+nocc 1911 nocc=nocc1*nocc2*nocc3*nocc4 1912 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%noccmmp,(nocc1,nocc2,nocc3,nocc4)) 1913 paw_ij_gathered(iat)%noccmmp(1:nocc1,1:nocc2,1:nocc3,1:nocc4)= & 1914 & reshape(buf_dp_all(indx_dp:indx_dp+nocc-1),(/nocc1,nocc2,nocc3,nocc4/)) 1915 indx_dp=indx_dp+nocc 1916 end if 1917 if (paw_ij_gathered(iat)%has_exexch_pot >=1) then 1918 sz1=buf_int_all(indx_int);indx_int=indx_int+1 1919 sz2=buf_int_all(indx_int);indx_int=indx_int+1 1920 sz3=buf_int_all(indx_int);indx_int=indx_int+1 1921 LIBPAW_ALLOCATE(paw_ij_gathered(iat)%vpawx,(sz1,sz2,sz3)) 1922 if (paw_ij_gathered(iat)%has_exexch_pot == 2) then 1923 paw_ij_gathered(iat)%vpawx(:,:,:)=& 1924 & reshape(buf_dp_all(indx_dp:indx_dp+sz1*sz2*sz3-1),(/sz1,sz2,sz3/)) 1925 indx_dp=indx_dp+sz1*sz2*sz3 1926 end if 1927 end if 1928 end do 1929 indx_int=indx_int-1;indx_dp=indx_dp-1 1930 if ((indx_int/=buf_int_size_all).or.(indx_dp/=buf_dp_size_all)) then 1931 write(msg,*) 'Wrong buffer sizes: buf_int_size_all=',buf_int_size_all, & 1932 & ' indx_int=',indx_int, ' buf_dp_size_all=',buf_dp_size_all,' indx_dp=',indx_dp 1933 LIBPAW_BUG(msg) 1934 end if 1935 end if 1936 1937 !Free memory 1938 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 1939 LIBPAW_DEALLOCATE(buf_int) 1940 LIBPAW_DEALLOCATE(buf_dp) 1941 LIBPAW_DEALLOCATE(buf_int_all) 1942 LIBPAW_DEALLOCATE(buf_dp_all) 1943 1944 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.
SOURCE
322 subroutine paw_ij_init(Paw_ij,cplex,nspinor,nsppol,nspden,pawspnorb,natom,ntypat,typat,Pawtab,& 323 & has_dij,has_dij0,has_dijfock,has_dijfr,has_dijhartree,has_dijhat,& ! Optional 324 & has_dijxc,has_dijxc_hat,has_dijxc_val,has_dijnd,has_dijso,has_dijU,has_dijexxc,& ! Optional 325 & has_exexch_pot,has_pawu_occ,nucdipmom,& ! Optional 326 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 327 328 !Arguments ------------------------------------ 329 !scalars 330 integer,intent(in) :: cplex,nspinor,nspden,nsppol,natom,ntypat,pawspnorb 331 integer,optional,intent(in) :: has_dij,has_dij0,has_dijfr,has_dijhat,has_dijxc,has_dijxc_hat,has_dijxc_val 332 integer,optional,intent(in) :: has_dijnd,has_dijso,has_dijhartree,has_dijfock,has_dijU,has_dijexxc 333 integer,optional,intent(in) :: has_exexch_pot,has_pawu_occ 334 integer,optional,intent(in) :: comm_atom 335 336 !arrays 337 integer,intent(in) :: typat(natom) 338 integer,optional,target,intent(in) :: mpi_atmtab(:) 339 real(dp),optional,intent(in) :: nucdipmom(3,natom) 340 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 341 type(Pawtab_type),intent(in) :: Pawtab(ntypat) 342 343 !Local variables------------------------------- 344 !scalars 345 integer :: cplex_dij,iat,iat_tot,itypat,lmn2_size,my_comm_atom,my_natom,ndij,qphase 346 logical :: my_atmtab_allocated,paral_atom,with_nucdipmom 347 !arrays 348 integer,pointer :: my_atmtab(:) 349 350 ! ************************************************************************* 351 352 !@Paw_ij_type 353 354 with_nucdipmom=.false.;if (present(nucdipmom)) with_nucdipmom=any(abs(nucdipmom)>tol8) 355 356 !Set up parallelism over atoms 357 my_natom=size(paw_ij);if (my_natom==0) return 358 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 359 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 360 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 361 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 362 363 do iat=1,my_natom 364 iat_tot=iat;if (paral_atom) iat_tot=my_atmtab(iat) 365 itypat=typat(iat_tot) 366 367 cplex_dij=1;if ((nspinor==2).or.(with_nucdipmom)) cplex_dij=2 368 qphase=cplex 369 370 lmn2_size =Pawtab(itypat)%lmn2_size 371 Paw_ij(iat)%qphase =qphase 372 Paw_ij(iat)%cplex_dij =cplex_dij 373 Paw_ij(iat)%itypat =itypat 374 Paw_ij(iat)%nspden =nspden 375 Paw_ij(iat)%nsppol =nsppol 376 Paw_ij(iat)%lmn_size =Pawtab(itypat)%lmn_size 377 Paw_ij(iat)%lmn2_size =lmn2_size 378 Paw_ij(iat)%ndij =MAX(nspinor**2,nspden) 379 380 ndij=Paw_ij(iat)%ndij 381 382 ! ================================== 383 ! === Allocations (all optional) === 384 ! ================================== 385 386 ! === Allocation for total Dij === 387 Paw_ij(iat)%has_dij=0 388 if (PRESENT(has_dij)) then 389 if (has_dij/=0) then 390 Paw_ij(iat)%has_dij=1 391 LIBPAW_ALLOCATE(Paw_ij(iat)%dij,(cplex_dij*qphase*lmn2_size,ndij)) 392 Paw_ij(iat)%dij(:,:)=zero 393 end if 394 end if 395 396 ! === Allocation for atomic Dij === 397 Paw_ij(iat)%has_dij0=0 398 if (PRESENT(has_dij0)) then 399 if (has_dij0/=0) then 400 Paw_ij(iat)%has_dij0=1 401 LIBPAW_ALLOCATE(Paw_ij(iat)%dij0,(lmn2_size)) 402 Paw_ij(iat)%dij0(:)=zero 403 end if 404 end if 405 406 ! === Allocation for Dij local exact exchange === 407 Paw_ij(iat)%has_dijexxc=0 408 if (PRESENT(has_dijexxc)) then 409 if (has_dijexxc/=0.and.Pawtab(itypat)%useexexch/=0) then 410 Paw_ij(iat)%has_dijexxc=1 411 LIBPAW_ALLOCATE(Paw_ij(iat)%dijexxc,(cplex_dij*lmn2_size,ndij)) 412 Paw_ij(iat)%dijexxc(:,:)=zero 413 end if 414 end if 415 416 ! === Allocation for Dij_Fock === 417 Paw_ij(iat)%has_dijfock=0 418 if (PRESENT(has_dijfock)) then 419 if (has_dijfock/=0) then 420 Paw_ij(iat)%has_dijfock=1 421 LIBPAW_ALLOCATE(Paw_ij(iat)%dijfock,(cplex_dij*lmn2_size,ndij)) 422 Paw_ij(iat)%dijfock(:,:)=zero 423 end if 424 end if 425 426 ! === Allocation for frozen part of 1st-order Dij === 427 Paw_ij(iat)%has_dijfr=0 428 if (PRESENT(has_dijfr)) then 429 if (has_dijfr/=0) then 430 Paw_ij(iat)%has_dijfr=1 431 LIBPAW_ALLOCATE(Paw_ij(iat)%dijfr,(cplex_dij*qphase*lmn2_size,ndij)) 432 Paw_ij(iat)%dijfr(:,:)=zero 433 end if 434 end if 435 436 ! === Allocation for Dij_Hartree === 437 Paw_ij(iat)%has_dijhartree=0 438 if (PRESENT(has_dijhartree)) then 439 if (has_dijhartree/=0) then 440 Paw_ij(iat)%has_dijhartree=1 441 LIBPAW_ALLOCATE(Paw_ij(iat)%dijhartree,(qphase*lmn2_size)) 442 Paw_ij(iat)%dijhartree(:)=zero 443 end if 444 end if 445 446 ! === Allocation for Dij_hat === 447 Paw_ij(iat)%has_dijhat=0 448 if (PRESENT(has_dijhat)) then 449 if (has_dijhat/=0) then 450 Paw_ij(iat)%has_dijhat=1 451 LIBPAW_ALLOCATE(Paw_ij(iat)%dijhat,(cplex_dij*qphase*lmn2_size,ndij)) 452 Paw_ij(iat)%dijhat(:,:)=zero 453 end if 454 end if 455 456 ! === Allocation for Dij nuclear dipole moment === 457 Paw_ij(iat)%has_dijnd=0 458 if (PRESENT(has_dijnd)) then 459 if (has_dijnd/=0.and.with_nucdipmom) then 460 Paw_ij(iat)%has_dijnd=1 461 LIBPAW_ALLOCATE(Paw_ij(iat)%dijnd,(cplex_dij*lmn2_size,ndij)) 462 Paw_ij(iat)%dijnd(:,:)=zero 463 end if 464 end if 465 466 ! === Allocation for Dij_SO === 467 Paw_ij(iat)%has_dijso=0 468 if (PRESENT(has_dijso)) then 469 if (has_dijso/=0.and.pawspnorb>0) then 470 Paw_ij(iat)%has_dijso=1 471 LIBPAW_ALLOCATE(Paw_ij(iat)%dijso,(cplex_dij*qphase*lmn2_size,ndij)) 472 Paw_ij(iat)%dijso(:,:)=zero 473 end if 474 end if 475 476 ! === Allocation for Dij_U_val === 477 Paw_ij(iat)%has_dijU=0 478 if (PRESENT(has_dijU)) then 479 if (has_dijU/=0) then 480 Paw_ij(iat)%has_dijU=1 481 LIBPAW_ALLOCATE(Paw_ij(iat)%dijU,(cplex_dij*qphase*lmn2_size,ndij)) 482 Paw_ij(iat)%dijU(:,:)=zero 483 end if 484 end if 485 486 ! === Allocation for total Dij_XC === 487 Paw_ij(iat)%has_dijxc=0 488 if (PRESENT(has_dijxc)) then 489 if (has_dijxc/=0) then 490 Paw_ij(iat)%has_dijxc=1 491 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc,(cplex_dij*qphase*lmn2_size,ndij)) 492 Paw_ij(iat)%dijxc(:,:)=zero 493 end if 494 end if 495 496 ! === Allocation for total Dij_XC_hat === 497 Paw_ij(iat)%has_dijxc_hat=0 498 if (PRESENT(has_dijxc_hat)) then 499 if (has_dijxc_hat/=0) then 500 Paw_ij(iat)%has_dijxc_hat=1 501 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc_hat,(cplex_dij*lmn2_size,ndij)) 502 Paw_ij(iat)%dijxc_hat(:,:)=zero 503 end if 504 end if 505 506 ! === Allocation for total Dij_XC_val === 507 Paw_ij(iat)%has_dijxc_val=0 508 if (PRESENT(has_dijxc_val)) then 509 if (has_dijxc_val/=0) then 510 Paw_ij(iat)%has_dijxc_val=1 511 LIBPAW_ALLOCATE(Paw_ij(iat)%dijxc_val,(cplex_dij*lmn2_size,ndij)) 512 Paw_ij(iat)%dijxc_val(:,:)=zero 513 end if 514 end if 515 516 ! === Allocation for PAW+U occupations === 517 Paw_ij(iat)%has_pawu_occ=0 518 if (PRESENT(has_pawu_occ)) then 519 if (has_pawu_occ/=0.and.Pawtab(itypat)%usepawu/=0) then 520 Paw_ij(iat)%has_pawu_occ=1 521 LIBPAW_ALLOCATE(Paw_ij(iat)%noccmmp,(cplex_dij,2*Pawtab(itypat)%lpawu+1,2*Pawtab(itypat)%lpawu+1,ndij)) 522 LIBPAW_ALLOCATE(Paw_ij(iat)%nocctot,(ndij)) 523 Paw_ij(iat)%noccmmp(:,:,:,:)=zero 524 Paw_ij(iat)%nocctot(:)=zero 525 end if 526 end if 527 528 ! === Allocation for PAW+LEXX potential === 529 Paw_ij(iat)%has_exexch_pot=0 530 if (PRESENT(has_exexch_pot)) then 531 if (has_exexch_pot/=0.and.Pawtab(itypat)%useexexch/=0) then 532 Paw_ij(iat)%has_exexch_pot=1 533 ! TODO solve issue with first dimension 534 LIBPAW_ALLOCATE(Paw_ij(iat)%vpawx,(1,lmn2_size,nspden)) 535 Paw_ij(iat)%vpawx(:,:,:)=zero 536 end if 537 end if 538 539 end do 540 541 !Destroy atom table used for parallelism 542 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 543 544 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
SOURCE
2400 subroutine paw_ij_isendreceive_getbuffer(paw_ij,npaw_ij_send,atm_indx_recv,buf_int,buf_dp) 2401 2402 !Arguments ------------------------------------ 2403 !scalars 2404 integer,intent(in) :: npaw_ij_send 2405 !arrays 2406 integer,intent(in):: atm_indx_recv(:),buf_int(:) 2407 real(dp),intent(in):: buf_dp(:) 2408 type(paw_ij_type),target,intent(inout) :: paw_ij(:) 2409 2410 !Local variables------------------------------- 2411 !scalars 2412 integer :: buf_dp_size,buf_int_size 2413 integer :: cplxq_lmn2_size,cplxdij_lmn2_size,cplxdijq_lmn2_size 2414 integer :: iat,iatom_tot,ii,ij,indx_dp,indx_int,ndij,nocc,nocc1,nocc2,nocc3,nocc4 2415 integer :: lmn2_size,sz1,sz2,sz3 2416 character(len=500) :: msg 2417 type(Paw_ij_type),pointer :: paw_ij1 2418 !arrays 2419 2420 ! ********************************************************************* 2421 2422 buf_int_size=size(buf_int) 2423 buf_dp_size=size(buf_dp) 2424 indx_int=1;indx_dp=1 2425 2426 do ij=1,npaw_ij_send 2427 iatom_tot=buf_int(indx_int) ;indx_int=indx_int+1 2428 iat= atm_indx_recv(iatom_tot) 2429 paw_ij1=>paw_ij(iat) 2430 paw_ij1%cplex_dij=buf_int(indx_int) ;indx_int=indx_int+1 2431 paw_ij1%qphase=buf_int(indx_int) ;indx_int=indx_int+1 2432 paw_ij1%itypat=buf_int(indx_int) ;indx_int=indx_int+1 2433 paw_ij1%nspden=buf_int(indx_int) ;indx_int=indx_int+1 2434 paw_ij1%nsppol=buf_int(indx_int) ;indx_int=indx_int+1 2435 paw_ij1%lmn_size=buf_int(indx_int) ;indx_int=indx_int+1 2436 paw_ij1%lmn2_size=buf_int(indx_int) ;indx_int=indx_int+1 2437 paw_ij1%ndij=buf_int(indx_int) ;indx_int=indx_int+1 2438 paw_ij1%has_dij=buf_int(indx_int) ;indx_int=indx_int+1 2439 paw_ij1%has_dij0=buf_int(indx_int) ;indx_int=indx_int+1 2440 paw_ij1%has_dijexxc=buf_int(indx_int) ;indx_int=indx_int+1 2441 paw_ij1%has_dijfock=buf_int(indx_int) ;indx_int=indx_int+1 2442 paw_ij1%has_dijfr=buf_int(indx_int) ;indx_int=indx_int+1 2443 paw_ij1%has_dijhartree=buf_int(indx_int) ;indx_int=indx_int+1 2444 paw_ij1%has_dijhat=buf_int(indx_int) ;indx_int=indx_int+1 2445 paw_ij1%has_dijnd=buf_int(indx_int) ;indx_int=indx_int+1 2446 paw_ij1%has_dijso=buf_int(indx_int) ;indx_int=indx_int+1 2447 paw_ij1%has_dijU=buf_int(indx_int) ;indx_int=indx_int+1 2448 paw_ij1%has_dijxc=buf_int(indx_int) ;indx_int=indx_int+1 2449 paw_ij1%has_dijxc_hat=buf_int(indx_int) ;indx_int=indx_int+1 2450 paw_ij1%has_dijxc_val=buf_int(indx_int) ;indx_int=indx_int+1 2451 paw_ij1%has_exexch_pot=buf_int(indx_int) ;indx_int=indx_int+1 2452 paw_ij1%has_pawu_occ=buf_int(indx_int) ;indx_int=indx_int+1 2453 if (paw_ij1%has_pawu_occ>=1) then 2454 nocc1=buf_int(indx_int) ;indx_int=indx_int+1 2455 nocc2=buf_int(indx_int) ;indx_int=indx_int+1 2456 nocc3=buf_int(indx_int) ;indx_int=indx_int+1 2457 nocc4=buf_int(indx_int) ;indx_int=indx_int+1 2458 else 2459 nocc1=0;nocc2=0;nocc3=0;nocc4=0 2460 end if 2461 lmn2_size=paw_ij1%lmn2_size 2462 cplxdij_lmn2_size=paw_ij1%cplex_dij*lmn2_size 2463 cplxq_lmn2_size=paw_ij1%qphase*lmn2_size 2464 cplxdijq_lmn2_size=cplxq_lmn2_size*paw_ij1%qphase 2465 ndij=paw_ij1%ndij 2466 2467 if (paw_ij1%has_dij>=1) then 2468 ii=cplxdijq_lmn2_size 2469 LIBPAW_ALLOCATE(paw_ij1%dij,(ii,ndij)) 2470 if (paw_ij1%has_dij==2) then 2471 paw_ij1%dij(:,:)= & 2472 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2473 indx_dp=indx_dp+ii*ndij 2474 end if 2475 end if 2476 if (paw_ij1%has_dij0 >=1) then 2477 ii=lmn2_size 2478 LIBPAW_ALLOCATE(paw_ij1%dij0,(ii)) 2479 if (paw_ij1%has_dij0==2) then 2480 paw_ij1%dij0(:)=buf_dp(indx_dp:indx_dp+ii-1) 2481 indx_dp=indx_dp+ii 2482 end if 2483 end if 2484 if (paw_ij1%has_dijexxc >=1) then 2485 ii=cplxdij_lmn2_size 2486 LIBPAW_ALLOCATE(paw_ij1%dijexxc,(ii,ndij)) 2487 if (paw_ij1%has_dijexxc==2) then 2488 paw_ij1%dijexxc(:,:)= & 2489 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2490 indx_dp=indx_dp+ii*ndij 2491 end if 2492 end if 2493 if (paw_ij1%has_dijfock >=1) then 2494 ii=cplxdij_lmn2_size 2495 LIBPAW_ALLOCATE(paw_ij1%dijfock,(ii,ndij)) 2496 if (paw_ij1%has_dijfock==2) then 2497 paw_ij1%dijfock(:,:)= & 2498 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2499 indx_dp=indx_dp+ii*ndij 2500 end if 2501 end if 2502 if (paw_ij1%has_dijfr >=1) then 2503 ii=cplxdijq_lmn2_size 2504 LIBPAW_ALLOCATE(paw_ij1%dijfr,(ii,ndij)) 2505 if (paw_ij1%has_dijfr==2) then 2506 paw_ij1%dijfr(:,:)= & 2507 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2508 indx_dp=indx_dp+ii*ndij 2509 end if 2510 end if 2511 if (paw_ij1%has_dijhartree >=1) then 2512 ii=cplxq_lmn2_size 2513 LIBPAW_ALLOCATE(paw_ij1%dijhartree,(ii)) 2514 if (paw_ij1%has_dijhartree==2) then 2515 paw_ij1%dijhartree(:)=buf_dp(indx_dp:indx_dp+ii-1) 2516 indx_dp=indx_dp+ii 2517 end if 2518 end if 2519 if (paw_ij1%has_dijhat >=1) then 2520 ii=cplxdijq_lmn2_size 2521 LIBPAW_ALLOCATE(paw_ij1%dijhat,(ii,ndij)) 2522 if (paw_ij1%has_dijhat==2) then 2523 paw_ij1%dijhat(:,:)= & 2524 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2525 indx_dp=indx_dp+ii*ndij 2526 end if 2527 end if 2528 if (paw_ij1%has_dijnd >=1) then 2529 ii=cplxdij_lmn2_size 2530 LIBPAW_ALLOCATE(paw_ij1%dijnd,(ii,ndij)) 2531 if (paw_ij1%has_dijnd==2) then 2532 paw_ij1%dijnd(:,:)= & 2533 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2534 indx_dp=indx_dp+ii*ndij 2535 end if 2536 end if 2537 if (paw_ij1%has_dijso >=1) then 2538 ii=cplxdijq_lmn2_size 2539 LIBPAW_ALLOCATE(paw_ij1%dijso,(ii,ndij)) 2540 if (paw_ij1%has_dijso==2) then 2541 paw_ij1%dijso(:,:)= & 2542 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2543 indx_dp=indx_dp+ii*ndij 2544 end if 2545 end if 2546 if (paw_ij1%has_dijU >=1) then 2547 ii=cplxdijq_lmn2_size 2548 LIBPAW_ALLOCATE(paw_ij1%dijU,(ii,ndij)) 2549 if (paw_ij1%has_dijU==2) then 2550 paw_ij1%dijU(:,:)= & 2551 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2552 indx_dp=indx_dp+ii*ndij 2553 end if 2554 end if 2555 if (paw_ij1%has_dijxc >=1) then 2556 ii=cplxdijq_lmn2_size 2557 LIBPAW_ALLOCATE(paw_ij1%dijxc,(ii,ndij)) 2558 if (paw_ij1%has_dijxc==2) then 2559 paw_ij1%dijxc(:,:)= & 2560 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2561 indx_dp=indx_dp+ii*ndij 2562 end if 2563 end if 2564 if (paw_ij1%has_dijxc_hat >=1) then 2565 ii=cplxdij_lmn2_size 2566 LIBPAW_ALLOCATE(paw_ij1%dijxc_hat,(ii,ndij)) 2567 if (paw_ij1%has_dijxc_hat==2) then 2568 paw_ij1%dijxc_hat(:,:)= & 2569 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2570 indx_dp=indx_dp+ii*ndij 2571 end if 2572 end if 2573 if (paw_ij1%has_dijxc_val >=1) then 2574 ii=cplxdij_lmn2_size 2575 LIBPAW_ALLOCATE(paw_ij1%dijxc_val,(ii,ndij)) 2576 if (paw_ij1%has_dijxc_val==2) then 2577 paw_ij1%dijxc_val(:,:)= & 2578 & reshape(buf_dp(indx_dp:indx_dp+ii*ndij-1),(/ii,ndij/)) 2579 indx_dp=indx_dp+ii*ndij 2580 end if 2581 end if 2582 if (paw_ij1%has_pawu_occ >=1) then 2583 nocc=paw_ij1%ndij 2584 LIBPAW_ALLOCATE(paw_ij1%nocctot,(nocc)) 2585 paw_ij1%nocctot(1:nocc)=buf_dp(indx_dp:indx_dp+nocc-1) 2586 indx_dp=indx_dp+nocc 2587 nocc=nocc1*nocc2*nocc3*nocc4 2588 LIBPAW_ALLOCATE(paw_ij1%noccmmp,(nocc1,nocc2,nocc3,nocc4)) 2589 paw_ij1%noccmmp(1:nocc1,1:nocc2,1:nocc3,1:nocc4)= & 2590 & reshape(buf_dp(indx_dp:indx_dp+nocc-1),(/nocc1,nocc2,nocc3,nocc4/)) 2591 indx_dp=indx_dp+nocc 2592 end if 2593 if (paw_ij1%has_exexch_pot >=1) then 2594 sz1=buf_int(indx_int);indx_int=indx_int+1 2595 sz2=buf_int(indx_int);indx_int=indx_int+1 2596 sz3=buf_int(indx_int);indx_int=indx_int+1 2597 LIBPAW_ALLOCATE(paw_ij1%vpawx,(sz1,sz2,sz3)) 2598 if (paw_ij1%has_exexch_pot == 2) then 2599 paw_ij1%vpawx(:,:,:)=& 2600 & reshape(buf_dp(indx_dp:indx_dp+sz1*sz2*sz3-1),(/sz1,sz2,sz3/)) 2601 indx_dp=indx_dp+sz1*sz2*sz3 2602 end if 2603 end if 2604 end do 2605 2606 if ((indx_int/=1+buf_int_size).or.(indx_dp /=1+buf_dp_size)) then 2607 write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size 2608 LIBPAW_BUG(msg) 2609 end if 2610 2611 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.
SOURCE
663 subroutine paw_ij_nullify(Paw_ij) 664 665 !Arguments ------------------------------------ 666 !arrays 667 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 668 669 !Local variables------------------------------- 670 integer :: iat,natom 671 672 ! ************************************************************************* 673 674 !@Paw_ij_type 675 676 ! MGPAW: This one could be removed/renamed, 677 ! variables can be initialized in the datatype declaration 678 ! Do we need to expose this in the public API? 679 680 natom=SIZE(Paw_ij(:));if (natom==0) return 681 682 ! Set all has_* flags to zero. 683 do iat=1,natom 684 ! === Set all has_* flags to zero === 685 Paw_ij(iat)%has_dij =0 686 Paw_ij(iat)%has_dij0 =0 687 Paw_ij(iat)%has_dijexxc =0 688 Paw_ij(iat)%has_dijfock =0 689 Paw_ij(iat)%has_dijfr =0 690 Paw_ij(iat)%has_dijhartree=0 691 Paw_ij(iat)%has_dijhat =0 692 Paw_ij(iat)%has_dijnd =0 693 Paw_ij(iat)%has_dijso =0 694 Paw_ij(iat)%has_dijU =0 695 Paw_ij(iat)%has_dijxc =0 696 Paw_ij(iat)%has_dijxc_hat =0 697 Paw_ij(iat)%has_dijxc_val =0 698 Paw_ij(iat)%has_exexch_pot=0 699 Paw_ij(iat)%has_pawu_occ =0 700 end do !iat 701 702 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
SOURCE
958 subroutine paw_ij_print(Paw_ij,unit,pawprtvol,pawspnorb,mode_paral,enunit,ipert, & 959 & mpi_atmtab,comm_atom,natom) 960 961 !Arguments ------------------------------------ 962 !scalars 963 integer,optional,intent(in) :: enunit,ipert 964 integer,optional,intent(in) :: comm_atom,natom 965 integer,optional,intent(in) :: pawprtvol,pawspnorb 966 integer,optional,intent(in) :: unit 967 character(len=4),optional,intent(in) :: mode_paral 968 !arrays 969 integer,optional,target,intent(in) :: mpi_atmtab(:) 970 type(Paw_ij_type),target,intent(in) :: Paw_ij(:) 971 972 !Local variables------------------------------- 973 character(len=7),parameter :: dspin(6)=(/"up ","down ","up-up ","dwn-dwn","up-dwn ","dwn-up "/) 974 !scalars 975 integer :: cplex_dij,iatom,iatom_tot,idij,idij_sym,lmn2_size,lmn_size,my_comm_atom,my_natom,nspden !klmn, 976 integer :: nsploop,nsppol,my_unt,ndij,qphase,tmp_cplex_dij,my_ipert,my_enunit,my_prtvol,size_paw_ij 977 logical :: my_atmtab_allocated,paral_atom 978 character(len=4) :: my_mode 979 character(len=2000) :: msg 980 !arrays 981 integer :: idum(0) 982 integer,pointer :: my_atmtab(:) 983 real(dp),allocatable,target :: dij(:),dijs(:),dijh(:,:) 984 real(dp),pointer :: dij2p(:),dij2p_(:) 985 986 ! ************************************************************************* 987 988 if (.False.) write(std_out,*)"pawspnorb:",pawspnorb 989 990 !@Paw_ij_type 991 size_paw_ij=SIZE(Paw_ij);if (size_paw_ij==0) return 992 993 my_unt =std_out ; if (PRESENT(unit )) my_unt =unit 994 my_prtvol=0 ; if (PRESENT(pawprtvol )) my_prtvol=pawprtvol 995 my_mode ='COLL' ; if (PRESENT(mode_paral)) my_mode =mode_paral 996 my_ipert =0 ; if (PRESENT(ipert)) my_ipert =ipert 997 my_enunit=0 ; if (PRESENT(enunit)) my_enunit=enunit 998 my_natom=size_paw_ij; if (PRESENT(natom)) my_natom=natom 999 1000 !Set up parallelism over atoms 1001 paral_atom=(present(comm_atom).and.my_natom/=size_paw_ij) 1002 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1003 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 1004 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,my_natom,my_natom_ref=size_paw_ij) 1005 1006 if (abs(my_prtvol)>=1) then 1007 if (my_ipert==0) then 1008 write(msg,'(4a)')ch10,' ==== Values of psp strength Dij (Hartree) ============' 1009 else 1010 write(msg,'(4a)')ch10,' ==== Values of psp strength Dij(1) (Hartree) =========' 1011 end if 1012 call wrtout(my_unt,msg,my_mode) 1013 end if 1014 1015 nsppol = Paw_ij(1)%nsppol 1016 nspden = Paw_ij(1)%nspden 1017 nsploop= nsppol; if (Paw_ij(1)%ndij==4) nsploop=4 1018 1019 do iatom=1,size_paw_ij 1020 1021 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1022 1023 lmn_size = Paw_ij(iatom)%lmn_size 1024 lmn2_size = Paw_ij(iatom)%lmn2_size 1025 cplex_dij = Paw_ij(iatom)%cplex_dij 1026 qphase = Paw_ij(iatom)%qphase 1027 ndij = Paw_ij(iatom)%ndij 1028 1029 ! ==================================== 1030 ! === Loop over density components === 1031 ! ==================================== 1032 do idij=1,nsploop 1033 1034 idij_sym=idij;if (ndij==4.and.idij>2) idij_sym=7-idij 1035 if (qphase==2) then 1036 LIBPAW_ALLOCATE(dij,(2*lmn2_size)) 1037 LIBPAW_ALLOCATE(dijs,(2*lmn2_size)) 1038 end if 1039 1040 ! =================== Detailed output ===================================== 1041 if (ABS(my_prtvol)>=1.and.(iatom_tot==1.or.iatom_tot==my_natom.or.my_prtvol<0)) then 1042 1043 !Title 1044 if (nspden==2.and.nsppol==1) then 1045 write(msg,'(2a,i3,3a)')ch10,& 1046 & ' >>>>>>>>>> Atom ',iatom_tot,':',ch10,& 1047 & ' (antiferromagnetism case: only one spin component)' 1048 else if (paw_ij(iatom)%ndij==1) then 1049 write(msg, '(2a,i3,a)') ch10,& 1050 & ' >>>>>>>>>> Atom ',iatom_tot,':' 1051 else 1052 write(msg,'(2a,i3,3a)') ch10,& 1053 & ' >>>>>>>>>> Atom ',iatom_tot,' (component ',TRIM(dspin(idij+2*(nsploop/4))),'):' 1054 end if 1055 call wrtout(my_unt,msg,my_mode) 1056 1057 !Dij atomic 1058 if (Paw_ij(iatom)%has_dij0/=0.and.idij<=2.and.my_ipert<=0) then 1059 write(msg,'(a)') ' ************ Dij atomic (Dij0) ***********' 1060 call wrtout(my_unt,msg,my_mode) 1061 call pawio_print_ij(my_unt,Paw_ij(iatom)%dij0,lmn2_size,1,lmn_size,-1,idum,0,my_prtvol,idum,-1.d0,1,& 1062 & opt_sym=2,mode_paral=my_mode) 1063 end if 1064 1065 !Dij Local Exact Exchange 1066 if (Paw_ij(iatom)%has_dijexxc/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1067 write(msg,'(a)') ' ************* Dij_Local Exact exchange **********' 1068 call wrtout(my_unt,msg,my_mode) 1069 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijexxc) 1070 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1071 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1072 end if 1073 1074 !Dij Fock 1075 if (Paw_ij(iatom)%has_dijfock/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1076 write(msg,'(a)') ' ************* Dij_Fock **********' 1077 call wrtout(my_unt,msg,my_mode) 1078 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijfock) 1079 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1080 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1081 end if 1082 1083 !Dij Frozen (RF) 1084 if (Paw_ij(iatom)%has_dijfr/=0.and.(idij<=2.or.nspden==4).and.my_ipert>0) then 1085 write(msg,'(a)') ' ************** Dij(1) Frozen **************' 1086 call wrtout(my_unt,msg,my_mode) 1087 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dijfr) 1088 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1089 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1090 end if 1091 1092 !Dij Hartree 1093 if (Paw_ij(iatom)%has_dijhartree/=0.and.idij<=2) then 1094 if (my_ipert==0) then 1095 write(msg,'(a)') ' ************** Dij Hartree ***************' 1096 else 1097 write(msg,'(a)') ' ************* Dij(1) Hartree *************' 1098 end if 1099 call wrtout(my_unt,msg,my_mode) 1100 LIBPAW_ALLOCATE(dijh,(qphase*lmn2_size,1)) 1101 dijh(:,1)=Paw_ij(iatom)%dijhartree(:) 1102 call get_dij_parts(1,qphase,dijh) 1103 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0, & 1104 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1105 LIBPAW_DEALLOCATE(dijh) 1106 end if 1107 1108 !Dij Hat 1109 if (Paw_ij(iatom)%has_dijhat/=0.and.(idij<=2.or.nspden==4)) then 1110 if (my_ipert==0) then 1111 write(msg,'(a)') ' **************** Dij_hat *****************' 1112 else 1113 write(msg,'(a)') ' ***** Dij_hat(1) (incl. frozen Dij) ******' 1114 end if 1115 call wrtout(my_unt,msg,my_mode) 1116 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dijhat) 1117 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1118 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1119 end if 1120 1121 !Dij nuclear dipole 1122 if (Paw_ij(iatom)%has_dijnd/=0) then 1123 write(msg,'(a)') ' *********** Dij Nuclear Dipole **********' 1124 call wrtout(my_unt,msg,my_mode) 1125 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijnd,always_img=.true.) 1126 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1127 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1128 end if 1129 1130 !Dij spin-orbit 1131 if (Paw_ij(iatom)%has_dijso/=0.and.my_ipert<=0) then 1132 write(msg,'(a)') ' ************** Dij SpinOrbit ************' 1133 call wrtout(my_unt,msg,my_mode) 1134 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dijso,always_img=.true.) 1135 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1136 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1137 end if 1138 1139 !Dij DFT+U 1140 if (Paw_ij(iatom)%has_dijU/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1141 write(msg,'(a)') ' ************* Dij_DFT+U (dijpawu) **********' 1142 call wrtout(my_unt,msg,my_mode) 1143 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%diju) 1144 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1145 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1146 end if 1147 1148 !Dij XC 1149 if (Paw_ij(iatom)%has_dijxc/=0.and.(idij<=2.or.nspden==4)) then 1150 if (my_ipert<=0) then 1151 write(msg,'(a)') ' ***************** Dij_xc *****************' 1152 else 1153 write(msg,'(a)') ' **************** Dij(1)_xc ***************' 1154 end if 1155 call wrtout(my_unt,msg,my_mode) 1156 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dijxc) 1157 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1158 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1159 end if 1160 1161 !Dij hat XC 1162 if (Paw_ij(iatom)%has_dijxc_hat/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1163 if (my_ipert<=0) then 1164 write(msg,'(a)') ' *************** Dijhat_xc ****************' 1165 else 1166 write(msg,'(a)') ' ************** Dij(1)hat_xc **************' 1167 end if 1168 call wrtout(my_unt,msg,my_mode) 1169 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijxc_hat) 1170 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1171 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1172 end if 1173 1174 !Dij XC val 1175 if (Paw_ij(iatom)%has_dijxc_val/=0.and.(idij<=2.or.nspden==4).and.my_ipert<=0) then 1176 write(msg, '(a)') ' *************** Dij_xc_val ***************' 1177 call wrtout(my_unt,msg,my_mode) 1178 call get_dij_parts(cplex_dij,1,Paw_ij(iatom)%dijxc_val) 1179 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1180 & my_prtvol,idum,-1.d0,1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1181 end if 1182 1183 !Dij TOTAL 1184 if (Paw_ij(iatom)%has_dij/=0) then 1185 if (my_ipert<=0) then 1186 write(msg,'(a)') ' ********** TOTAL Dij in Ha **********' 1187 else 1188 write(msg,'(a)') ' ********** TOTAL Dij(1) in Ha **********' 1189 end if 1190 call wrtout(my_unt,msg,my_mode) 1191 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dij,always_img=.true.) 1192 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1193 & my_prtvol,idum,50.d0*dble(3-2*idij),1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1194 if (my_enunit>0) then 1195 if (my_ipert<=0) then 1196 write(msg,'(a)') ' ********** TOTAL Dij in eV **********' 1197 else 1198 write(msg,'(a)') ' ********** TOTAL Dij(1) in eV **********' 1199 end if 1200 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1201 & my_prtvol,idum,-1._dp,2,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1202 end if 1203 end if 1204 1205 if (allocated(dij)) then 1206 LIBPAW_DEALLOCATE(dij) 1207 end if 1208 if (allocated(dijs)) then 1209 LIBPAW_DEALLOCATE(dijs) 1210 end if 1211 1212 end if !(ABS(my_prtvol)>=1.and.(iatom_tot==1.or.iatom_tot==my_natom.or.my_prtvol<0) 1213 1214 ! =================== Standard output ===================================== 1215 if ((abs(my_prtvol)==0).and.(iatom_tot==1.or.iatom_tot==my_natom)) then 1216 1217 !Title 1218 if (idij==1) then 1219 if (my_ipert<=0) then 1220 write(msg, '(2a,i6,a)') ch10,' ****** Psp strength Dij in Ha (atom ',iatom_tot,') *****' 1221 else 1222 write(msg, '(2a,i6,a)') ch10,' **** Psp strength Dij(1) in Ha (atom ',iatom_tot,') *****' 1223 end if 1224 if (nspden==2.and.nsppol==1) then 1225 write(msg,'(4a)') trim(msg),') *****',ch10,' (antiferromagnetism case: only one spin component)' 1226 end if 1227 call wrtout(my_unt,msg,my_mode) 1228 end if 1229 if (paw_ij(iatom)%ndij/=1) then 1230 write(msg,'(3a)') ' Component ',trim(dspin(idij+2*(nsploop/4))),':' 1231 call wrtout(my_unt,msg,my_mode) 1232 end if 1233 1234 !Dij TOTAL 1235 if (Paw_ij(iatom)%has_dij/=0) then 1236 call get_dij_parts(cplex_dij,qphase,Paw_ij(iatom)%dij,always_img=.true.) 1237 call pawio_print_ij(my_unt,dij2p,lmn2_size,tmp_cplex_dij,lmn_size,-1,idum,0,& 1238 & my_prtvol,idum,50.d0*dble(3-2*idij),1,opt_sym=2,asym_ij=dij2p_,mode_paral=my_mode) 1239 end if 1240 1241 end if 1242 1243 ! =================== End main loops ===================================== 1244 1245 if (allocated(dij)) then 1246 LIBPAW_DEALLOCATE(dij) 1247 end if 1248 if (allocated(dijs)) then 1249 LIBPAW_DEALLOCATE(dijs) 1250 end if 1251 1252 end do !idij 1253 end do !iat 1254 1255 call wrtout(my_unt,' ',my_mode) 1256 1257 !Small helper function 1258 contains 1259 1260 !Real and imaginary parts of phase. 1261 subroutine get_dij_parts(my_cplex_dij,my_qphase,my_dij,always_img) 1262 1263 integer,intent(in) :: my_cplex_dij,my_qphase 1264 logical,intent(in),optional :: always_img 1265 real(dp),intent(in),target :: my_dij(:,:) 1266 integer :: my_idij,my_idij_sym,kk 1267 logical :: always_img_ 1268 always_img_=.false.;if(present(always_img)) always_img_=always_img 1269 my_idij=min(size(my_dij,2),idij) 1270 my_idij_sym=min(size(my_dij,2),idij_sym) 1271 if (my_qphase==1) then 1272 if ((idij<=nsppol.or.idij==2).and.(.not.always_img_))then 1273 tmp_cplex_dij=1 1274 dij2p => my_dij(1:my_cplex_dij*lmn2_size:my_cplex_dij,my_idij) 1275 dij2p_ => dij2p 1276 else 1277 tmp_cplex_dij=my_cplex_dij 1278 dij2p => my_dij(1:my_cplex_dij*lmn2_size:1,my_idij) 1279 dij2p_ => my_dij(1:my_cplex_dij*lmn2_size:1,my_idij_sym) 1280 end if 1281 else 1282 tmp_cplex_dij=2 1283 if (my_cplex_dij==1) then 1284 do kk=1,lmn2_size 1285 dij(2*kk-1)= my_dij(kk,my_idij) 1286 dij(2*kk )= my_dij(kk+lmn2_size,my_idij) 1287 dijs(2*kk-1)= my_dij(kk,my_idij) 1288 dijs(2*kk )=-my_dij(kk+lmn2_size,my_idij) 1289 end do 1290 else 1291 do kk=1,lmn2_size 1292 dij(2*kk-1)= my_dij(2*kk-1,idij)-my_dij(2*kk +2*lmn2_size,my_idij) 1293 dij(2*kk )= my_dij(2*kk ,idij)+my_dij(2*kk-1+2*lmn2_size,my_idij) 1294 dijs(2*kk-1)= my_dij(2*kk-1,idij_sym)+my_dij(2*kk +2*lmn2_size,my_idij_sym) 1295 dijs(2*kk )= my_dij(2*kk ,idij_sym)-my_dij(2*kk-1+2*lmn2_size,my_idij_sym) 1296 end do 1297 end if 1298 dij2p => dij ; dij2p_ => dijs 1299 end if 1300 end subroutine get_dij_parts 1301 1302 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
SOURCE
1985 subroutine paw_ij_redistribute(paw_ij,mpi_comm_in,mpi_comm_out,& 1986 & natom,mpi_atmtab_in,mpi_atmtab_out,paw_ij_out,& 1987 & SendAtomProc,SendAtomList,RecvAtomProc,RecvAtomList) 1988 1989 !Arguments ------------------------------------ 1990 !scalars 1991 integer,intent(in) :: mpi_comm_in,mpi_comm_out 1992 integer,optional,intent(in) :: natom 1993 !arrays 1994 integer,intent(in),optional,target :: mpi_atmtab_in(:),mpi_atmtab_out(:) 1995 type(paw_ij_type),allocatable,intent(inout) :: paw_ij(:) 1996 type(paw_ij_type),pointer,optional :: paw_ij_out(:) !vz_i 1997 integer,intent(in),optional :: SendAtomProc(:),SendAtomList(:),RecvAtomProc(:),RecvAtomList(:) 1998 1999 !Local variables------------------------------- 2000 !scalars 2001 integer :: algo_option,i1,iatom,iat_in,iat_out,ierr,iircv,iisend,imsg,imsg_current 2002 integer :: imsg1,iproc_rcv,iproc_send,ireq,me_exch,mpi_comm_exch,my_natom_in,my_natom_out,my_tag,natom_tot 2003 integer :: nb_dp,nb_int,nb_msg,nbmsg_incoming,nbrecv,nbrecvmsg,nbsendreq,nbsent,nbsend,next,npaw_ij_sent 2004 integer :: nproc_in,nproc_out 2005 logical :: flag,in_place,message_yet_prepared,my_atmtab_in_allocated,my_atmtab_out_allocated,paral_atom 2006 !arrays 2007 integer :: buf_size(3),request1(3) 2008 integer,pointer :: my_atmtab_in(:),my_atmtab_out(:) 2009 integer,allocatable :: atmtab_send(:),atm_indx_in(:),atm_indx_out(:),From(:),buf_int1(:),request(:) 2010 integer,allocatable,target:: buf_int(:) 2011 integer,pointer :: buf_ints(:) 2012 logical,allocatable :: msg_pick(:) 2013 real(dp),allocatable :: buf_dp1(:) 2014 real(dp),allocatable,target :: buf_dp(:) 2015 real(dp),pointer :: buf_dps(:) 2016 type(coeffi1_type),target,allocatable :: tab_buf_int(:),tab_buf_atom(:) 2017 type(coeff1_type),target,allocatable :: tab_buf_dp(:) 2018 type(paw_ij_type),allocatable :: paw_ij_all(:) 2019 type(paw_ij_type),pointer :: paw_ij_out1(:) 2020 2021 ! ************************************************************************* 2022 2023 !@Paw_ij_type 2024 2025 in_place=(.not.present(paw_ij_out)) 2026 my_natom_in=size(paw_ij) 2027 2028 !If not "in_place", destroy the output datastructure 2029 if (.not.in_place) then 2030 if (associated(paw_ij_out)) then 2031 call paw_ij_free(paw_ij_out) 2032 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_out) 2033 end if 2034 end if 2035 2036 !Special sequential case 2037 if (mpi_comm_in==xmpi_comm_self.and.mpi_comm_out==xmpi_comm_self) then 2038 if ((.not.in_place).and.(my_natom_in>0)) then 2039 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_in)) 2040 call paw_ij_nullify(paw_ij_out) 2041 call paw_ij_copy(paw_ij,paw_ij_out) 2042 end if 2043 return 2044 end if 2045 2046 !Get total natom 2047 if (present(natom)) then 2048 natom_tot=natom 2049 else 2050 natom_tot=my_natom_in 2051 call xmpi_sum(natom_tot,mpi_comm_in,ierr) 2052 end if 2053 2054 !Select input distribution 2055 if (present(mpi_atmtab_in)) then 2056 my_atmtab_in => mpi_atmtab_in 2057 my_atmtab_in_allocated=.false. 2058 else 2059 call get_my_atmtab(mpi_comm_in,my_atmtab_in,my_atmtab_in_allocated,& 2060 & paral_atom,natom_tot,my_natom_in) 2061 end if 2062 2063 !Select output distribution 2064 if (present(mpi_atmtab_out)) then 2065 my_natom_out=size(mpi_atmtab_out) 2066 my_atmtab_out => mpi_atmtab_out 2067 my_atmtab_out_allocated=.false. 2068 else 2069 call get_my_natom(mpi_comm_out,my_natom_out,natom_tot) 2070 call get_my_atmtab(mpi_comm_out,my_atmtab_out,my_atmtab_out_allocated,& 2071 & paral_atom,natom_tot) 2072 end if 2073 2074 !Select algo according to optional input arguments 2075 algo_option=1 2076 if (present(SendAtomProc).and.present(SendAtomList).and.& 2077 & present(RecvAtomProc).and.present(RecvAtomList)) algo_option=2 2078 2079 2080 !Brute force algorithm (allgather + scatter) 2081 !--------------------------------------------------------- 2082 if (algo_option==1) then 2083 2084 LIBPAW_DATATYPE_ALLOCATE(paw_ij_all,(natom_tot)) 2085 call paw_ij_nullify(paw_ij_all) 2086 call paw_ij_copy(paw_ij,paw_ij_all,comm_atom=mpi_comm_in,mpi_atmtab=my_atmtab_in) 2087 if (in_place) then 2088 call paw_ij_free(paw_ij) 2089 LIBPAW_DATATYPE_DEALLOCATE(paw_ij) 2090 LIBPAW_DATATYPE_ALLOCATE(paw_ij,(my_natom_out)) 2091 call paw_ij_nullify(paw_ij) 2092 call paw_ij_copy(paw_ij_all,paw_ij,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 2093 else 2094 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_out)) 2095 call paw_ij_nullify(paw_ij_out) 2096 call paw_ij_copy(paw_ij_all,paw_ij_out,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 2097 end if 2098 call paw_ij_free(paw_ij_all) 2099 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_all) 2100 2101 2102 !Asynchronous algorithm (asynchronous communications) 2103 !--------------------------------------------------------- 2104 else if (algo_option==2) then 2105 2106 nbsend=size(SendAtomProc) ; nbrecv=size(RecvAtomProc) 2107 2108 if (in_place) then 2109 if (my_natom_out > 0) then 2110 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out1,(my_natom_out)) 2111 call paw_ij_nullify(paw_ij_out1) 2112 else 2113 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out1,(0)) 2114 end if 2115 else 2116 LIBPAW_DATATYPE_ALLOCATE(paw_ij_out,(my_natom_out)) 2117 call paw_ij_nullify(paw_ij_out) 2118 paw_ij_out1=>paw_ij_out 2119 end if 2120 2121 nproc_in=xmpi_comm_size(mpi_comm_in) 2122 nproc_out=xmpi_comm_size(mpi_comm_out) 2123 if (nproc_in<=nproc_out) mpi_comm_exch=mpi_comm_out 2124 if (nproc_in>nproc_out) mpi_comm_exch=mpi_comm_in 2125 me_exch=xmpi_comm_rank(mpi_comm_exch) 2126 2127 2128 ! Dimension put to the maximum to send 2129 LIBPAW_ALLOCATE(atmtab_send,(nbsend)) 2130 LIBPAW_ALLOCATE(atm_indx_in,(natom_tot)) 2131 atm_indx_in=-1 2132 do iatom=1,my_natom_in 2133 atm_indx_in(my_atmtab_in(iatom))=iatom 2134 end do 2135 LIBPAW_ALLOCATE(atm_indx_out,(natom_tot)) 2136 atm_indx_out=-1 2137 do iatom=1,my_natom_out 2138 atm_indx_out(my_atmtab_out(iatom))=iatom 2139 end do 2140 2141 LIBPAW_DATATYPE_ALLOCATE(tab_buf_int,(nbsend)) 2142 LIBPAW_DATATYPE_ALLOCATE(tab_buf_dp,(nbsend)) 2143 LIBPAW_DATATYPE_ALLOCATE(tab_buf_atom,(nbsend)) 2144 LIBPAW_ALLOCATE(request,(3*nbsend)) 2145 2146 ! A send buffer in an asynchrone communication couldn't be deallocate before it has been receive 2147 nbsent=0 ; ireq=0 ; iisend=0 ; nbsendreq=0 ; nb_msg=0 2148 do iisend=1,nbsend 2149 iproc_rcv=SendAtomProc(iisend) 2150 next=-1 2151 if (iisend < nbsend) next=SendAtomProc(iisend+1) 2152 if (iproc_rcv /= me_exch) then 2153 nbsent=nbsent+1 2154 atmtab_send(nbsent)=SendAtomList(iisend) ! we groups the atoms sends to the same process 2155 if (iproc_rcv /= next) then 2156 if (nbsent > 0 ) then 2157 ! Check if message has been yet prepared 2158 message_yet_prepared=.false. 2159 do imsg=1,nb_msg 2160 if (size(tab_buf_atom(imsg)%value) /= nbsent) then 2161 cycle 2162 else 2163 do imsg1=1,nbsent 2164 if (tab_buf_atom(imsg)%value(imsg1)/=atmtab_send(imsg1)) exit 2165 message_yet_prepared=.true. 2166 imsg_current=imsg 2167 end do 2168 end if 2169 enddo 2170 ! Create the message 2171 if (.not.message_yet_prepared) then 2172 nb_msg=nb_msg+1 2173 call paw_ij_isendreceive_fillbuffer( & 2174 & paw_ij,atmtab_send,atm_indx_in,nbsent,buf_int,nb_int,buf_dp,nb_dp) 2175 LIBPAW_ALLOCATE(tab_buf_int(nb_msg)%value,(nb_int)) 2176 LIBPAW_ALLOCATE(tab_buf_dp(nb_msg)%value,(nb_dp)) 2177 tab_buf_int(nb_msg)%value(1:nb_int)=buf_int(1:nb_int) 2178 tab_buf_dp(nb_msg)%value(1:nb_dp)=buf_dp(1:nb_dp) 2179 LIBPAW_DEALLOCATE(buf_int) 2180 LIBPAW_DEALLOCATE(buf_dp) 2181 LIBPAW_ALLOCATE(tab_buf_atom(nb_msg)%value, (nbsent)) 2182 tab_buf_atom(nb_msg)%value(1:nbsent)=atmtab_send(1:nbsent) 2183 imsg_current=nb_msg 2184 end if 2185 ! Communicate 2186 buf_size(1)=size(tab_buf_int(imsg_current)%value) 2187 buf_size(2)=size(tab_buf_dp(imsg_current)%value) 2188 buf_size(3)=nbsent 2189 buf_ints=>tab_buf_int(imsg_current)%value 2190 buf_dps=>tab_buf_dp(imsg_current)%value 2191 my_tag=200 2192 ireq=ireq+1 2193 call xmpi_isend(buf_size,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2194 my_tag=201 2195 ireq=ireq+1 2196 call xmpi_isend(buf_ints,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2197 my_tag=202 2198 ireq=ireq+1 2199 call xmpi_isend(buf_dps,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 2200 nbsendreq=ireq 2201 nbsent=0 2202 end if 2203 end if 2204 else ! Just a renumbering, not a sending 2205 iat_in=atm_indx_in(SendAtomList(iisend)) 2206 iat_out=atm_indx_out(my_atmtab_in(iat_in)) 2207 call paw_ij_copy(paw_ij(iat_in:iat_in),paw_ij_out1(iat_out:iat_out)) 2208 nbsent=0 2209 end if 2210 end do 2211 2212 LIBPAW_ALLOCATE(From,(nbrecv)) 2213 From(:)=-1 ; nbrecvmsg=0 2214 do iircv=1,nbrecv 2215 iproc_send=RecvAtomProc(iircv) !receive from (RcvAtomProc is sorted by growing process) 2216 next=-1 2217 if (iircv < nbrecv) next=RecvAtomProc(iircv+1) 2218 if (iproc_send /= me_exch .and. iproc_send/=next) then 2219 nbrecvmsg=nbrecvmsg+1 2220 From(nbrecvmsg)=iproc_send 2221 end if 2222 end do 2223 2224 LIBPAW_ALLOCATE(msg_pick,(nbrecvmsg)) 2225 msg_pick=.false. 2226 nbmsg_incoming=nbrecvmsg 2227 do while (nbmsg_incoming > 0) 2228 do i1=1,nbrecvmsg 2229 if (.not.msg_pick(i1)) then 2230 iproc_send=From(i1) 2231 flag=.false. 2232 my_tag=200 2233 call xmpi_iprobe(iproc_send,my_tag,mpi_comm_exch,flag,ierr) 2234 if (flag) then 2235 msg_pick(i1)=.true. 2236 call xmpi_irecv(buf_size,iproc_send,my_tag,mpi_comm_exch,request1(1),ierr) 2237 call xmpi_wait(request1(1),ierr) 2238 nb_int=buf_size(1) 2239 nb_dp=buf_size(2) 2240 npaw_ij_sent=buf_size(3) 2241 LIBPAW_ALLOCATE(buf_int1,(nb_int)) 2242 LIBPAW_ALLOCATE(buf_dp1,(nb_dp)) 2243 my_tag=201 2244 call xmpi_irecv(buf_int1,iproc_send,my_tag,mpi_comm_exch,request1(2),ierr) 2245 my_tag=202 2246 call xmpi_irecv(buf_dp1,iproc_send,my_tag,mpi_comm_exch,request1(3),ierr) 2247 call xmpi_waitall(request1(2:3),ierr) 2248 call paw_ij_isendreceive_getbuffer(paw_ij_out1,npaw_ij_sent,atm_indx_out,buf_int1,buf_dp1) 2249 nbmsg_incoming=nbmsg_incoming-1 2250 LIBPAW_DEALLOCATE(buf_int1) 2251 LIBPAW_DEALLOCATE(buf_dp1) 2252 end if 2253 end if 2254 end do 2255 end do 2256 LIBPAW_DEALLOCATE(msg_pick) 2257 2258 if (in_place) then 2259 call paw_ij_free(paw_ij) 2260 LIBPAW_DATATYPE_DEALLOCATE(paw_ij) 2261 LIBPAW_DATATYPE_ALLOCATE(paw_ij,(my_natom_out)) 2262 call paw_ij_copy(paw_ij_out1,paw_ij) 2263 call paw_ij_free(paw_ij_out1) 2264 LIBPAW_DATATYPE_DEALLOCATE(paw_ij_out1) 2265 end if 2266 2267 ! Wait for deallocating arrays that all sending operations has been realized 2268 if (nbsendreq > 0) then 2269 call xmpi_waitall(request(1:nbsendreq),ierr) 2270 end if 2271 2272 ! Deallocate buffers 2273 do i1=1,nb_msg 2274 LIBPAW_DEALLOCATE(tab_buf_int(i1)%value) 2275 LIBPAW_DEALLOCATE(tab_buf_dp(i1)%value) 2276 LIBPAW_DEALLOCATE(tab_buf_atom(i1)%value) 2277 end do 2278 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_int) 2279 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_dp) 2280 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_atom) 2281 LIBPAW_DEALLOCATE(From) 2282 LIBPAW_DEALLOCATE(request) 2283 LIBPAW_DEALLOCATE(atmtab_send) 2284 LIBPAW_DEALLOCATE(atm_indx_in) 2285 LIBPAW_DEALLOCATE(atm_indx_out) 2286 2287 end if !algo_option 2288 2289 !Eventually release temporary pointers 2290 call free_my_atmtab(my_atmtab_in,my_atmtab_in_allocated) 2291 call free_my_atmtab(my_atmtab_out,my_atmtab_out_allocated) 2292 2293 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
SOURCE
2310 subroutine paw_ij_reset_flags(Paw_ij,all,dijhartree,self_consistent) 2311 2312 !Arguments ------------------------------------ 2313 !scalars 2314 logical,optional,intent(in) :: all,dijhartree,self_consistent 2315 !arrays 2316 type(Paw_ij_type),intent(inout) :: Paw_ij(:) 2317 2318 !Local variables------------------------------- 2319 integer :: iat,natom 2320 logical :: all_,dijhartree_,self_consistent_ 2321 2322 ! ************************************************************************* 2323 2324 !@Paw_ij_type 2325 2326 natom=SIZE(Paw_ij);if (natom==0) return 2327 all_=.true.;if (present(all)) all_=all 2328 dijhartree_=.false.;if (present(dijhartree)) dijhartree_=dijhartree 2329 self_consistent_=.false.;if (present(self_consistent)) self_consistent_=self_consistent 2330 2331 if (dijhartree_) then 2332 do iat=1,natom 2333 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2334 end do 2335 2336 else if (self_consistent_) then 2337 do iat=1,natom 2338 if (Paw_ij(iat)%has_dij >0) Paw_ij(iat)%has_dij =1 2339 if (Paw_ij(iat)%has_dijexxc >0) Paw_ij(iat)%has_dijexxc =1 2340 if (Paw_ij(iat)%has_dijfock >0) Paw_ij(iat)%has_dijfock =1 2341 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2342 if (Paw_ij(iat)%has_dijhat >0) Paw_ij(iat)%has_dijhat =1 2343 if (Paw_ij(iat)%has_dijnd >0) Paw_ij(iat)%has_dijnd =1 2344 if (Paw_ij(iat)%has_dijso >0) Paw_ij(iat)%has_dijso =1 2345 if (Paw_ij(iat)%has_dijU >0) Paw_ij(iat)%has_dijU =1 2346 if (Paw_ij(iat)%has_dijxc >0) Paw_ij(iat)%has_dijxc =1 2347 if (Paw_ij(iat)%has_dijxc_hat >0) Paw_ij(iat)%has_dijxc_hat =1 2348 if (Paw_ij(iat)%has_dijxc_val >0) Paw_ij(iat)%has_dijxc_val =1 2349 if (Paw_ij(iat)%has_exexch_pot>0) Paw_ij(iat)%has_exexch_pot=1 2350 if (Paw_ij(iat)%has_pawu_occ >0) Paw_ij(iat)%has_pawu_occ =1 2351 end do 2352 2353 else if (all_) then 2354 do iat=1,natom 2355 if (Paw_ij(iat)%has_dij >0) Paw_ij(iat)%has_dij =1 2356 if (Paw_ij(iat)%has_dij0 >0) Paw_ij(iat)%has_dij0 =1 2357 if (Paw_ij(iat)%has_dijexxc >0) Paw_ij(iat)%has_dijexxc =1 2358 if (Paw_ij(iat)%has_dijfock >0) Paw_ij(iat)%has_dijfock =1 2359 if (Paw_ij(iat)%has_dijfr >0) Paw_ij(iat)%has_dijfr =1 2360 if (Paw_ij(iat)%has_dijhartree>0) Paw_ij(iat)%has_dijhartree=1 2361 if (Paw_ij(iat)%has_dijhat >0) Paw_ij(iat)%has_dijhat =1 2362 if (Paw_ij(iat)%has_dijnd >0) Paw_ij(iat)%has_dijnd =1 2363 if (Paw_ij(iat)%has_dijso >0) Paw_ij(iat)%has_dijso =1 2364 if (Paw_ij(iat)%has_dijU >0) Paw_ij(iat)%has_dijU =1 2365 if (Paw_ij(iat)%has_dijxc >0) Paw_ij(iat)%has_dijxc =1 2366 if (Paw_ij(iat)%has_dijxc_hat >0) Paw_ij(iat)%has_dijxc_hat =1 2367 if (Paw_ij(iat)%has_dijxc_val >0) Paw_ij(iat)%has_dijxc_val =1 2368 if (Paw_ij(iat)%has_exexch_pot>0) Paw_ij(iat)%has_exexch_pot=1 2369 if (Paw_ij(iat)%has_pawu_occ >0) Paw_ij(iat)%has_pawu_occ =1 2370 end do 2371 end if 2372 2373 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
52 type,public :: paw_ij_type 53 54 !Integer scalars 55 56 integer :: cplex_dij 57 ! cplex_dij=1 if dij are real 58 ! cplex_dij=2 if dij are complex (spin-orbit, non-collinear magnetism, magnetic field, ...) 59 60 integer :: has_dij=0 61 ! 1 if dij is allocated 62 ! 2 if dij is already computed 63 64 integer :: has_dij0=0 65 ! 1 if dij0 is allocated 66 ! 2 if dij0 is already computed 67 68 integer :: has_dijexxc=0 69 ! 1 if dijexxc is associated and used, 0 otherwise 70 ! 2 if dijexxc is already computed 71 72 integer :: has_dijfock=0 73 ! 1 if dijfock is allocated 74 ! 2 if dijfock is already computed 75 76 integer :: has_dijfr=0 77 ! 1 if dijfr is allocated 78 ! 2 if dijfr is already computed 79 80 integer :: has_dijhartree=0 81 ! 1 if dijhartree is allocated 82 ! 2 if dijhartree is already computed 83 84 integer :: has_dijhat=0 85 ! 1 if dijhat is allocated 86 ! 2 if dijhat is already computed 87 88 integer :: has_dijnd=0 89 ! on site term due to nuclear dipole moment 90 ! 1 if dijnd is associated and used, 0 otherwise 91 ! 2 if dijnd is already computed 92 93 integer :: has_dijso=0 94 ! 1 if dijso is associated and used, 0 otherwise 95 ! 2 if dijso is already computed 96 97 integer :: has_dijU=0 98 ! 1 if dijU is associated and used, 0 otherwise 99 ! 2 if dijU is already computed 100 101 integer :: has_dijxc=0 102 ! 1 if dijxc is associated and used, 0 otherwise 103 ! 2 if dijxc is already computed 104 105 integer :: has_dijxc_hat=0 106 ! 1 if dijxc_hat is associated and used, 0 otherwise 107 ! 2 if dijxc_hat is already computed 108 109 integer :: has_dijxc_val=0 110 ! 1 if dijxc_val is associated and used, 0 otherwise 111 ! 2 if dijxc_val is already computed 112 113 integer :: has_exexch_pot=0 114 ! 1 if PAW+(local exact exchange) potential is allocated 115 116 integer :: has_pawu_occ=0 117 ! 1 if PAW+U occupations are allocated 118 119 integer :: itypat 120 ! itypat=type of the atom 121 122 integer :: lmn_size 123 ! Number of (l,m,n) elements for the paw basis 124 125 integer :: lmn2_size 126 ! lmn2_size=lmn_size*(lmn_size+1)/2 127 ! where lmn_size is the number of (l,m,n) elements for the paw basis 128 129 integer :: ndij 130 ! Number of components of dij 131 ! Usually ndij=nspden, except for nspinor==2 (where ndij=nspinor**2) 132 133 integer :: nspden 134 ! Number of spin-density components (may be different from dtset%nspden if spin-orbit) 135 136 integer :: nsppol 137 ! Number of independant spin-components 138 139 integer :: qphase 140 ! qphase=2 if dij contain a exp(-i.q.r) phase (as in the q<>0 RF case), 1 if not 141 ! (this may change the ij symmetry) 142 143 !Real (real(dp)) arrays 144 145 real(dp), allocatable :: dij(:,:) 146 ! dij(cplex_dij*qphase*lmn2_size,ndij) 147 ! Dij term (non-local operator) 148 ! May be complex if cplex_dij=2 or qphase=2 149 ! ==== Storage for the 1st dimension ==== 150 ! For each klmn=ij: 151 ! When Dij is complex (cplex_dij=2): 152 ! dij(2*ij-1,:) contains the real part 153 ! dij(2*ij ,:) contains the imaginary part 154 ! When a exp(-i.q.r) phase is included (qphase=2): 155 ! dij(1:cplex_dij*lmn2_size,:) 156 ! contains the real part of the phase, i.e. D_ij*cos(q.r) 157 ! dij(cplex_dij*lmn2_size+1:2*cplex_dij*lmn2_size,:) 158 ! contains the imaginary part of the phase, i.e. D_ij*sin(q.r) 159 ! ==== Storage for the 2nd dimension ==== 160 ! dij(:,1) contains Dij^up-up 161 ! dij(:,2) contains Dij^dn-dn 162 ! dij(:,3) contains Dij^up-dn (only if nspinor=2) 163 ! dij(:,4) contains Dij^dn-up (only if nspinor=2) 164 165 real(dp), allocatable :: dij0(:) 166 ! dij0(lmn2_size) 167 ! Atomic part of Dij (read from PAW dataset) 168 ! Same storage as Dij (see above); always real, spin independent 169 170 real(dp), allocatable :: dijexxc(:,:) 171 ! dijexxc(cplex_dij*lmn2_size,ndij) 172 ! On-site matrix elements of the Fock operator (Local Exact exchange implementation) 173 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 174 175 real(dp), allocatable :: dijfock(:,:) 176 ! dijfock(cplex_dij*lmn2_size,ndij) 177 ! Dij_fock term 178 ! Contains all contributions to Dij from Fock exchange 179 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 180 181 real(dp), allocatable :: dijfr(:,:) 182 ! dijhat(cplex_dij*qphase*lmn2_size,ndij) 183 ! For response function calculation only 184 ! RF Frozen part of Dij (depends on q vector but not on 1st-order wave function) 185 ! Same storage as Dij (see above) 186 187 real(dp), allocatable :: dijhartree(:) 188 ! dijhartree(qphase*lmn2_size) 189 ! Dij_hartree term; contains all contributions to Dij from hartree 190 ! Warning: dimensioned only by qphase (exp(-iqr)), not cplex_dij 191 ! Same storage as Dij (see above); spin independent 192 193 real(dp), allocatable :: dijhat(:,:) 194 ! dijhat(cplex_dij*qphase*lmn2_size,ndij) 195 ! Dij_hat term (non-local operator) i.e \sum_LM \int_FFT Q_{ij}^{LM} vtrial 196 ! Same storage as Dij (see above) 197 ! Same storage as Dij (see above) 198 199 real(dp), allocatable :: dijnd(:,:) 200 ! dijnd(cplex_dij*lmn2_size,ndij) 201 ! On-site matrix elements of -\frac{1}{c}\mu\cdot L/r^3 202 ! Same storage as Dij (see above) 203 204 real(dp), allocatable :: dijso(:,:) 205 ! dijso(cplex_dij*qphase*lmn2_size,ndij) 206 ! On-site matrix elements of L.S i.e <phi_i|L.S|phi_j> 207 ! Same storage as Dij (see above) 208 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 209 210 real(dp), allocatable :: dijU(:,:) 211 ! dijU(cplex_dij*qphase*lmn2_size,ndij) 212 ! On-site matrix elements of the U part of the PAW Hamiltonian. 213 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 214 215 real(dp), allocatable :: dijxc(:,:) 216 ! dijxc(cplex_dij*qphase*lmn2_size,ndij) 217 ! On-site matrix elements of vxc i.e 218 ! <phi_i|vxc[n1+nc]|phi_j> - <tphi_i|vxc(tn1+nhat+tnc]|tphi_j> 219 ! Same storage as Dij (see above) 220 221 real(dp), allocatable :: dijxc_hat(:,:) 222 ! dijxc_hat(cplex_dij*lmn2_size,ndij) 223 ! Dij_hat term i.e \sum_LM \int_FFT Q_{ij}^{LM} Vxc 224 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 225 226 real(dp), allocatable :: dijxc_val(:,:) 227 ! dijxc_val(cplex_dij*lmn2_size,ndij) 228 ! Onsite matrix elements of valence-only vxc i.e 229 ! <phi_i|vxc[n1]|phi_j> - <tphi_i|vxc(tn1+nhat]|tphi_j> 230 ! Same storage as Dij (see above); not available for RF (i.e. qphase=2) 231 232 real(dp), allocatable :: noccmmp(:,:,:,:) 233 ! noccmmp(cplex_dij,2*lpawu+1,2*lpawu+1,nocc_nspden) 234 ! cplex_dij=1 if collinear 235 ! cplex_dij=2 if spin orbit is used 236 ! cplex_dij=2 is used if non-collinear (for coherence, it is not necessary in this case, however) 237 ! gives occupation matrix for lda+u (computed in setnoccmmp) 238 ! Stored as: noccmmp(:,:,1)= n^{up,up}_{m,mp} 239 ! noccmmp(:,:,2)= n^{dn,dn}_{m,mp} 240 ! noccmmp(:,:,3)= n^{up,dn}_{m,mp} 241 ! noccmmp(:,:,4)= n^{dn,up}_{m,mp} 242 ! noccmmp(m,mp,:) is computed from rhoij(klmn) with m=klmntomn(2)>mp=klmntomn(1) 243 244 real(dp), allocatable :: nocctot(:) 245 ! nocctot(ndij) 246 ! gives trace of occupation matrix for lda+u (computed in pawdenpot) 247 ! for each value of ispden (1 or 2) 248 249 real(dp), allocatable :: vpawx(:,:,:) 250 ! vpawx(1,2*lexexch+1,nspden) 251 ! exact exchange potential 252 253 end type paw_ij_type 254 255 !public procedures. 256 public :: paw_ij_init ! Creation method 257 public :: paw_ij_free ! Free memory 258 public :: paw_ij_nullify 259 public :: paw_ij_copy ! Copy object 260 public :: paw_ij_print ! Printout of the object 261 public :: paw_ij_gather ! MPI gather 262 public :: paw_ij_redistribute ! MPI redistribute 263 public :: paw_ij_reset_flags ! Resent the internal flags. 264 265 !private procedures. 266 private :: paw_ij_isendreceive_getbuffer 267 private :: paw_ij_isendreceive_fillbuffer