TABLE OF CONTENTS
- ABINIT/m_paw_an
- m_paw_an/paw_an_copy
- m_paw_an/paw_an_free
- m_paw_an/paw_an_gather
- m_paw_an/paw_an_init
- m_paw_an/paw_an_isendreceive_fillbuffer
- m_paw_an/paw_an_isendreceive_getbuffer
- m_paw_an/paw_an_nullify
- m_paw_an/paw_an_print
- m_paw_an/paw_an_redistribute
- m_paw_an/paw_an_reset_flags
- m_paw_an/paw_an_type
ABINIT/m_paw_an [ Modules ]
NAME
m_paw_an
FUNCTION
This module contains the definition of the paw_an_type structured datatype, as well as related functions and methods. paw_an_type variables contain various arrays given on ANgular mesh or ANgular moments for a given atom.
COPYRIGHT
Copyright (C) 2013-2018 ABINIT group (MT, FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
24 #include "libpaw.h" 25 26 MODULE m_paw_an 27 28 USE_DEFS 29 USE_MSG_HANDLING 30 USE_MPI_WRAPPERS 31 USE_MEMORY_PROFILING 32 33 use m_paral_atom, only : get_my_atmtab, free_my_atmtab, get_my_natom 34 use m_pawang, only : pawang_type 35 use m_pawtab, only : pawtab_type 36 37 implicit none 38 39 private 40 41 !public procedures. 42 public :: paw_an_init 43 public :: paw_an_free 44 public :: paw_an_nullify 45 public :: paw_an_copy 46 public :: paw_an_print 47 public :: paw_an_gather 48 public :: paw_an_redistribute 49 public :: paw_an_reset_flags 50 51 !private procedures. 52 private :: paw_an_isendreceive_getbuffer 53 private :: paw_an_isendreceive_fillbuffer
m_paw_an/paw_an_copy [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_copy
FUNCTION
Copy one paw_an datastructure into another Can take into accound changes of dimensions Can copy a shared paw_an 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_an_in(:)<type(paw_an_type)>= input paw_an datastructure
SIDE EFFECTS
paw_an_cpy(:)<type(paw_an_type)>= output paw_an datastructure
NOTES
paw_an_cpy must have been allocated in the calling function.
PARENTS
m_paw_an
CHILDREN
SOURCE
561 subroutine paw_an_copy(paw_an_in,paw_an_cpy,& 562 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 563 564 565 !This section has been created automatically by the script Abilint (TD). 566 !Do not modify the following lines by hand. 567 #undef ABI_FUNC 568 #define ABI_FUNC 'paw_an_copy' 569 !End of the abilint section 570 571 implicit none 572 573 !Arguments ------------------------------------ 574 !scalars 575 integer,optional,intent(in) :: comm_atom 576 !arrays 577 integer,optional,target,intent(in) :: mpi_atmtab(:) 578 type(Paw_an_type),intent(in),target :: paw_an_in(:) 579 type(Paw_an_type),intent(inout),target :: paw_an_cpy(:) 580 581 !Local variables------------------------------- 582 !scalars 583 integer :: cplx_mesh_size,ij,ij1,lm_size,my_comm_atom,my_natom,nkxc1,nk3xc1,npaw_an_in 584 integer :: npaw_an_max,npaw_an_out,nspden,paral_case,v_size,sz1 585 logical :: my_atmtab_allocated,paral_atom 586 character(len=500) :: msg 587 type(Paw_an_type),pointer :: paw_an_in1, paw_an_out1 588 !arrays 589 integer,pointer :: my_atmtab(:) 590 type(Paw_an_type), pointer :: paw_an_out(:) 591 592 ! ************************************************************************* 593 594 !@Paw_an_type 595 596 !Retrieve sizes 597 npaw_an_in=size(paw_an_in);npaw_an_out=size(paw_an_cpy) 598 599 !Set up parallelism over atoms 600 paral_atom=(present(comm_atom));if (paral_atom) paral_atom=(xmpi_comm_size(comm_atom)>1) 601 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 602 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 603 my_atmtab_allocated=.false. 604 605 !Determine in which case we are (parallelism, ...) 606 !No parallelism: a single copy operation 607 paral_case=0;npaw_an_max=npaw_an_in 608 paw_an_out => paw_an_cpy 609 if (paral_atom) then 610 if (npaw_an_out<npaw_an_in) then ! Parallelism: the copy operation is a scatter 611 call get_my_natom(my_comm_atom,my_natom,npaw_an_in) 612 if (my_natom==npaw_an_out) then 613 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,npaw_an_in) 614 paral_case=1;npaw_an_max=npaw_an_out 615 paw_an_out => paw_an_cpy 616 else 617 msg=' npaw_an_out should be equal to my_natom !' 618 MSG_BUG(msg) 619 end if 620 else ! Parallelism: the copy operation is a gather 621 call get_my_natom(my_comm_atom,my_natom,npaw_an_out) 622 if (my_natom==npaw_an_in) then 623 paral_case=2;npaw_an_max=npaw_an_in 624 else 625 msg=' npaw_ij_in should be equal to my_natom !' 626 MSG_BUG(msg) 627 end if 628 end if 629 end if 630 631 !First case: a simple copy or a scatter 632 if (npaw_an_max>0.and.((paral_case==0).or.(paral_case==1))) then 633 call paw_an_free(paw_an_cpy) 634 call paw_an_nullify(paw_an_cpy) 635 636 ! Loop on paw_ij components 637 do ij1=1,npaw_an_max 638 ij=ij1; if (paral_case==1) ij=my_atmtab(ij1) 639 640 paw_an_in1=>paw_an_in(ij) 641 paw_an_out1=>paw_an_out(ij1) 642 paw_an_out1%angl_size =paw_an_in1%angl_size 643 paw_an_out1%cplex =paw_an_in1%cplex 644 paw_an_out1%has_kxc =paw_an_in1%has_kxc 645 paw_an_out1%has_k3xc =paw_an_in1%has_k3xc 646 paw_an_out1%has_vhartree =paw_an_in1%has_vhartree 647 paw_an_out1%has_vxc =paw_an_in1%has_vxc 648 paw_an_out1%has_vxcval =paw_an_in1%has_vxcval 649 paw_an_out1%has_vxc_ex =paw_an_in1%has_vxc_ex 650 paw_an_out1%itypat =paw_an_in1%itypat 651 paw_an_out1%lm_size =paw_an_in1%lm_size 652 paw_an_out1%mesh_size =paw_an_in1%mesh_size 653 paw_an_out1%nkxc1 =paw_an_in1%nkxc1 654 paw_an_out1%nk3xc1 =paw_an_in1%nk3xc1 655 paw_an_out1%nspden =paw_an_in1%nspden 656 if (allocated(paw_an_in1%lmselect)) then 657 sz1=size(paw_an_in1%lmselect) 658 LIBPAW_ALLOCATE(paw_an_out1%lmselect,(sz1)) 659 paw_an_out1%lmselect(:)=paw_an_in1%lmselect(:) 660 end if 661 v_size=0 662 if (paw_an_in1%has_vxc>0) then 663 v_size=size(paw_an_in1%vxc1,2) 664 else if (paw_an_in1%has_kxc>0) then 665 v_size=size(paw_an_in1%kxc1,2) 666 else if (paw_an_in1%has_k3xc>0) then 667 v_size=size(paw_an_in1%k3xc1,2) 668 else if (paw_an_in1%has_vxcval>0) then 669 v_size=size(paw_an_in1%vxc1_val,2) 670 else if (paw_an_in1%has_vxc_ex>0) then 671 v_size=size(paw_an_in1%vxc_ex,2) 672 else if (paw_an_in1%has_vhartree>0) then 673 v_size=size(paw_an_in1%vh1,2) 674 end if 675 nspden=paw_an_in1%nspden 676 lm_size=paw_an_in1%lm_size 677 cplx_mesh_size=paw_an_in1%cplex*paw_an_in1%mesh_size 678 nkxc1=paw_an_in1%nkxc1 679 if (paw_an_in1%has_kxc>0) then 680 LIBPAW_ALLOCATE(paw_an_out1%kxc1,(cplx_mesh_size,v_size,nkxc1)) 681 LIBPAW_ALLOCATE(paw_an_out1%kxct1,(cplx_mesh_size,v_size,nkxc1)) 682 if (paw_an_in1%has_kxc==2.and.nkxc1>0) then 683 paw_an_out1%kxc1(:,:,:)=paw_an_in1%kxc1(:,:,:) 684 paw_an_out1%kxct1(:,:,:)=paw_an_in1%kxct1(:,:,:) 685 end if 686 end if 687 nk3xc1=paw_an_in1%nk3xc1 688 if (paw_an_in1%has_k3xc>0) then 689 LIBPAW_ALLOCATE(paw_an_out1%k3xc1,(cplx_mesh_size,v_size,nk3xc1)) 690 LIBPAW_ALLOCATE(paw_an_out1%k3xct1,(cplx_mesh_size,v_size,nk3xc1)) 691 if (paw_an_in1%has_k3xc==2.and.nk3xc1>0) then 692 paw_an_out1%k3xc1(:,:,:)=paw_an_in1%k3xc1(:,:,:) 693 paw_an_out1%k3xct1(:,:,:)=paw_an_in1%k3xct1(:,:,:) 694 end if 695 end if 696 if (paw_an_in1%has_vhartree>0) then 697 LIBPAW_ALLOCATE(paw_an_out1%vh1,(cplx_mesh_size,lm_size,nspden)) 698 LIBPAW_ALLOCATE(paw_an_out1%vht1,(cplx_mesh_size,lm_size,nspden)) 699 if (paw_an_in1%has_vhartree==2) then 700 paw_an_out1%vh1(:,:,:)=paw_an_in1%vh1(:,:,:) 701 paw_an_out1%vht1(:,:,:)=paw_an_in1%vht1(:,:,:) 702 end if 703 end if 704 if (paw_an_in1%has_vxc>0) then 705 LIBPAW_ALLOCATE(paw_an_out1%vxc1,(cplx_mesh_size,v_size,nspden)) 706 LIBPAW_ALLOCATE(paw_an_out1%vxct1,(cplx_mesh_size,v_size,nspden)) 707 if (paw_an_in1%has_vxc==2) then 708 paw_an_out1%vxc1(:,:,:)=paw_an_in1%vxc1(:,:,:) 709 paw_an_out1%vxct1(:,:,:)=paw_an_in1%vxct1(:,:,:) 710 end if 711 end if 712 if (paw_an_in1%has_vxcval>0) then 713 LIBPAW_ALLOCATE(paw_an_out1%vxc1_val,(cplx_mesh_size,v_size,nspden)) 714 LIBPAW_ALLOCATE(paw_an_out1%vxct1_val,(cplx_mesh_size,v_size,nspden)) 715 if (paw_an_in1%has_vxcval==2) then 716 paw_an_out1%vxc1_val(:,:,:)=paw_an_in1%vxc1_val(:,:,:) 717 paw_an_out1%vxct1_val(:,:,:)=paw_an_in1%vxct1_val(:,:,:) 718 end if 719 end if 720 if (paw_an_in1%has_vxc_ex>0) then 721 LIBPAW_ALLOCATE(paw_an_out1%vxc_ex,(cplx_mesh_size,v_size,nspden)) 722 if (paw_an_in1%has_vxc_ex==2) then 723 paw_an_out1%vxc_ex(:,:,:)=paw_an_in1%vxc_ex(:,:,:) 724 end if 725 end if 726 end do 727 end if 728 729 !Second case: a gather 730 if (paral_case==2) then 731 call paw_an_gather(paw_an_in,paw_an_cpy,-1,my_comm_atom,my_atmtab) 732 end if 733 734 !Destroy atom table used for parallelism 735 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 736 737 end subroutine paw_an_copy
m_paw_an/paw_an_free [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_free
FUNCTION
Deallocate pointers and nullify flags in a paw_an structure
SIDE EFFECTS
Paw_an(:)<type(Paw_an_type)>=various arrays given on ANgular mesh or ANgular moments All associated pointers in Paw_an(:) are deallocated
PARENTS
bethe_salpeter,dfpt_scfcv,m_paral_pert,m_paw_an,respfn,scfcv,screening sigma
CHILDREN
SOURCE
396 subroutine paw_an_free(Paw_an) 397 398 399 !This section has been created automatically by the script Abilint (TD). 400 !Do not modify the following lines by hand. 401 #undef ABI_FUNC 402 #define ABI_FUNC 'paw_an_free' 403 !End of the abilint section 404 405 implicit none 406 407 !Arguments ------------------------------------ 408 !arrays 409 type(Paw_an_type),intent(inout) :: Paw_an(:) 410 411 !Local variables------------------------------- 412 integer :: iat,natom 413 414 ! ************************************************************************* 415 416 !@Paw_an_type 417 418 natom=SIZE(Paw_an);if (natom==0) return 419 420 do iat=1,natom 421 if (allocated(Paw_an(iat)%lmselect )) then 422 LIBPAW_DEALLOCATE(Paw_an(iat)%lmselect) 423 end if 424 if (allocated(Paw_an(iat)%vh1 )) then 425 LIBPAW_DEALLOCATE(Paw_an(iat)%vh1) 426 end if 427 if (allocated(Paw_an(iat)%vht1 )) then 428 LIBPAW_DEALLOCATE(Paw_an(iat)%vht1) 429 end if 430 if (allocated(Paw_an(iat)%vxc1 )) then 431 LIBPAW_DEALLOCATE(Paw_an(iat)%vxc1) 432 end if 433 if (allocated(Paw_an(iat)%vxc1_val )) then 434 LIBPAW_DEALLOCATE(Paw_an(iat)%vxc1_val) 435 end if 436 if (allocated(Paw_an(iat)%vxct1 )) then 437 LIBPAW_DEALLOCATE(Paw_an(iat)%vxct1) 438 end if 439 if (allocated(Paw_an(iat)%vxct1_val)) then 440 LIBPAW_DEALLOCATE(Paw_an(iat)%vxct1_val) 441 end if 442 if (allocated(Paw_an(iat)%kxc1 )) then 443 LIBPAW_DEALLOCATE(Paw_an(iat)%kxc1) 444 end if 445 if (allocated(Paw_an(iat)%kxct1 )) then 446 LIBPAW_DEALLOCATE(Paw_an(iat)%kxct1) 447 end if 448 if (allocated(Paw_an(iat)%k3xc1 )) then 449 LIBPAW_DEALLOCATE(Paw_an(iat)%k3xc1) 450 end if 451 if (allocated(Paw_an(iat)%k3xct1 )) then 452 LIBPAW_DEALLOCATE(Paw_an(iat)%k3xct1) 453 end if 454 if (allocated(Paw_an(iat)%vxc_ex )) then 455 LIBPAW_DEALLOCATE(Paw_an(iat)%vxc_ex) 456 end if 457 458 ! === Reset all has_* flags === 459 Paw_an(iat)%has_kxc =0 460 Paw_an(iat)%has_k3xc =0 461 Paw_an(iat)%has_vhartree=0 462 Paw_an(iat)%has_vxc =0 463 Paw_an(iat)%has_vxcval =0 464 Paw_an(iat)%has_vxc_ex =0 465 end do !iat 466 467 end subroutine paw_an_free
m_paw_an/paw_an_gather [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_gather
FUNCTION
(All)Gather paw_an datastructures
INPUTS
master=master communicator receiving data ; if -1 do a ALLGATHER comm_atom= communicator over atom mpi_atmtab(:)=--optional-- indexes of the atoms treated by current calling proc paw_an_in(:)<type(paw_an_type)>= input paw_an datastructures on every process
OUTPUT
paw_an_gathered(:)<type(paw_an_type)>= output paw_an datastructure
PARENTS
m_paw_an
CHILDREN
SOURCE
878 subroutine paw_an_gather(Paw_an_in,paw_an_gathered,master,comm_atom,mpi_atmtab) 879 880 881 !This section has been created automatically by the script Abilint (TD). 882 !Do not modify the following lines by hand. 883 #undef ABI_FUNC 884 #define ABI_FUNC 'paw_an_gather' 885 !End of the abilint section 886 887 implicit none 888 889 !Arguments ------------------------------------ 890 integer,intent(in) :: master,comm_atom 891 !arrays 892 integer,optional,target,intent(in) :: mpi_atmtab(:) 893 type(Paw_an_type),target,intent(in) :: Paw_an_in(:) 894 type(Paw_an_type),target,intent(inout) :: Paw_an_gathered(:) 895 896 !Local variables------------------------------- 897 !scalars 898 integer :: buf_dp_size,buf_dp_size_all,buf_int_size,buf_int_size_all,cplx_mesh_size 899 integer :: iat,iatot,ierr,has_lm_select,i1,i2,ij,indx_int,indx_dp 900 integer :: lm_size,me_atom 901 integer :: my_natom,natom,nkxc1,nk3xc1,npaw_an_in_sum,nproc_atom,nspden,v_size,sz1,sz2,sz3 902 logical :: my_atmtab_allocated,paral_atom 903 character(len=500) :: msg 904 type(Paw_an_type),pointer :: paw_an_in1,paw_an_gathered1 905 !arrays 906 integer :: bufsz(2) 907 integer,allocatable :: buf_int(:),buf_int_all(:) 908 integer,allocatable :: count_dp(:),count_int(:),count_tot(:),displ_dp(:),displ_int(:) 909 integer,pointer :: my_atmtab(:) 910 real(dp),allocatable :: buf_dp(:),buf_dp_all(:) 911 912 ! ************************************************************************* 913 914 !@Paw_an_type 915 916 if (master/=-1) then 917 msg='simple gather (master/=-1) not yet implemented !' 918 MSG_BUG(msg) 919 end if 920 921 my_natom=size(paw_an_in);natom=size(paw_an_gathered) 922 923 !Set up parallelism over atoms 924 paral_atom=(my_natom/=natom) 925 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 926 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 927 nproc_atom=xmpi_comm_size(comm_atom) 928 me_atom=xmpi_comm_rank(comm_atom) 929 930 !Special case: one process (simple copy) 931 if (nproc_atom==1) then 932 if (master==-1.or.me_atom==master) then 933 call paw_an_free(paw_an_gathered) 934 call paw_an_nullify(paw_an_gathered) 935 do iat=1,my_natom 936 paw_an_in1=>paw_an_in(iat) 937 paw_an_gathered1%itypat =paw_an_in1%itypat 938 paw_an_gathered1%nspden =paw_an_in1%nspden 939 paw_an_gathered1%cplex =paw_an_in1%cplex 940 paw_an_gathered1%mesh_size =paw_an_in1%mesh_size 941 paw_an_gathered1%angl_size =paw_an_in1%angl_size 942 paw_an_gathered1%lm_size =paw_an_in1%lm_size 943 paw_an_gathered1%nkxc1 =paw_an_in1%nkxc1 944 paw_an_gathered1%nk3xc1 =paw_an_in1%nk3xc1 945 paw_an_gathered1%has_vxc =paw_an_in1%has_vxc 946 paw_an_gathered1%has_kxc =paw_an_in1%has_kxc 947 paw_an_gathered1%has_k3xc =paw_an_in1%has_k3xc 948 paw_an_gathered1%has_vxcval =paw_an_in1%has_vxcval 949 paw_an_gathered1%has_vxc_ex =paw_an_in1%has_vxc_ex 950 paw_an_gathered1%has_vhartree =paw_an_in1%has_vhartree 951 if (allocated(paw_an_in1%lmselect)) then 952 sz1=size(paw_an_in1%lmselect) 953 LIBPAW_ALLOCATE(paw_an_gathered1%lmselect,(sz1)) 954 paw_an_gathered1%lmselect(:)=paw_an_in1%lmselect(:) 955 end if 956 if (allocated(paw_an_in1%vxc1)) then 957 sz1=size(paw_an_in1%vxc1,1);sz2=size(paw_an_in1%vxc1,2) 958 sz3=size(paw_an_in1%vxc1,3) 959 LIBPAW_ALLOCATE(paw_an_gathered1%vxc1,(sz1,sz2,sz3)) 960 paw_an_gathered1%vxc1(:,:,:)=paw_an_in1%vxc1(:,:,:) 961 end if 962 if (allocated(paw_an_in1%vxct1)) then 963 sz1=size(paw_an_in1%vxct1,1);sz2=size(paw_an_in1%vxct1,2) 964 sz3=size(paw_an_in1%vxct1,3) 965 LIBPAW_ALLOCATE(paw_an_gathered1%vxct1,(sz1,sz2,sz3)) 966 paw_an_gathered1%vxct1(:,:,:)=paw_an_in1%vxct1(:,:,:) 967 end if 968 if (allocated(paw_an_in1%kxc1)) then 969 sz1=size(paw_an_in1%kxc1,1);sz2=size(paw_an_in1%kxc1,2) 970 sz3=size(paw_an_in1%kxc1,3) 971 LIBPAW_ALLOCATE(paw_an_gathered1%kxc1,(sz1,sz2,sz3)) 972 if (sz3>0) paw_an_gathered1%kxc1(:,:,:)=paw_an_in1%kxc1(:,:,:) 973 end if 974 if (allocated(paw_an_in1%k3xc1)) then 975 sz1=size(paw_an_in1%k3xc1,1);sz2=size(paw_an_in1%k3xc1,2) 976 sz3=size(paw_an_in1%k3xc1,3) 977 LIBPAW_ALLOCATE(paw_an_gathered1%k3xc1,(sz1,sz2,sz3)) 978 if (sz3>0) paw_an_gathered1%k3xc1(:,:,:)=paw_an_in1%k3xc1(:,:,:) 979 end if 980 if (allocated(paw_an_in1%kxct1)) then 981 sz1=size(paw_an_in1%kxct1,1);sz2=size(paw_an_in1%kxct1,2) 982 sz3=size(paw_an_in1%kxct1,3) 983 LIBPAW_ALLOCATE(paw_an_gathered1%kxct1,(sz1,sz2,sz3)) 984 if (sz3>0) paw_an_gathered1%kxct1(:,:,:)=paw_an_in1%kxct1(:,:,:) 985 end if 986 if (allocated(paw_an_in1%k3xct1)) then 987 sz1=size(paw_an_in1%k3xct1,1);sz2=size(paw_an_in1%k3xct1,2) 988 sz3=size(paw_an_in1%k3xct1,3) 989 LIBPAW_ALLOCATE(paw_an_gathered1%k3xct1,(sz1,sz2,sz3)) 990 if (sz3>0) paw_an_gathered1%k3xct1(:,:,:)=paw_an_in1%k3xct1(:,:,:) 991 end if 992 if (allocated(paw_an_in1%vxc1_val)) then 993 sz1=size(paw_an_in1%vxc1_val,1);sz2=size(paw_an_in1%vxc1_val,2) 994 sz3=size(paw_an_in1%vxc1_val,3) 995 LIBPAW_ALLOCATE(paw_an_gathered1%vxc1_val,(sz1,sz2,sz3)) 996 paw_an_gathered1%vxc1_val(:,:,:)=paw_an_in1%vxc1_val(:,:,:) 997 end if 998 if (allocated(paw_an_in1%vxct1_val)) then 999 sz1=size(paw_an_in1%vxct1_val,1);sz2=size(paw_an_in1%vxct1_val,2) 1000 sz3=size(paw_an_in1%vxct1_val,3) 1001 LIBPAW_ALLOCATE(paw_an_gathered1%vxct1_val,(sz1,sz2,sz3)) 1002 paw_an_gathered1%vxct1_val(:,:,:)=paw_an_in1%vxct1_val(:,:,:) 1003 end if 1004 if (allocated(paw_an_in1%vxc_ex)) then 1005 sz1=size(paw_an_in1%vxc_ex,1);sz2=size(paw_an_in1%vxc_ex,2) 1006 sz3=size(paw_an_in1%vxc_ex,3) 1007 LIBPAW_ALLOCATE(paw_an_gathered1%vxc_ex,(sz1,sz2,sz3)) 1008 paw_an_gathered1%vxc_ex(:,:,:)=paw_an_in1%vxc_ex(:,:,:) 1009 end if 1010 if (allocated(paw_an_in1%vh1)) then 1011 sz1=size(paw_an_in1%vh1,1);sz2=size(paw_an_in1%vh1,2) 1012 sz3=size(paw_an_in1%vh1,3) 1013 LIBPAW_ALLOCATE(paw_an_gathered1%vh1,(sz1,sz2,sz3)) 1014 paw_an_gathered1%vh1(:,:,:)=paw_an_in1%vh1(:,:,:) 1015 end if 1016 if (allocated(paw_an_in1%vht1)) then 1017 sz1=size(paw_an_in1%vht1,1);sz2=size(paw_an_in1%vht1,2) 1018 sz3=size(paw_an_in1%vht1,3) 1019 LIBPAW_ALLOCATE(paw_an_gathered1%vht1,(sz1,sz2,sz3)) 1020 paw_an_gathered1%vht1(:,:,:)=paw_an_in1%vht1(:,:,:) 1021 end if 1022 end do 1023 end if 1024 return 1025 end if 1026 1027 !Test on sizes 1028 npaw_an_in_sum=my_natom 1029 call xmpi_sum(npaw_an_in_sum,comm_atom,ierr) 1030 if (master==-1) then 1031 if (natom/=npaw_an_in_sum) then 1032 msg='Wrong sizes sum[npaw_an_in]/=natom !' 1033 MSG_BUG(msg) 1034 end if 1035 else 1036 if (me_atom==master.and.natom/=npaw_an_in_sum) then 1037 msg='(2) paw_an_gathered wrongly allocated !' 1038 MSG_BUG(msg) 1039 end if 1040 end if 1041 1042 !Compute sizes of buffers 1043 buf_int_size=0;buf_dp_size=0 1044 do ij=1,my_natom 1045 buf_int_size=buf_int_size+17+size(paw_an_in(ij)%lmselect) 1046 end do 1047 do ij=1,my_natom 1048 paw_an_in1=>paw_an_in(ij) 1049 if (paw_an_in1%has_vxc==2) then 1050 buf_dp_size=buf_dp_size+size(paw_an_in1%vxc1) 1051 buf_dp_size=buf_dp_size+size(paw_an_in1%vxct1) 1052 end if 1053 if (paw_an_in1%has_kxc==2) then 1054 buf_dp_size=buf_dp_size+size(paw_an_in1%kxc1) 1055 buf_dp_size=buf_dp_size+size(paw_an_in1%kxct1) 1056 end if 1057 if (paw_an_in1%has_k3xc==2) then 1058 buf_dp_size=buf_dp_size+size(paw_an_in1%k3xc1) 1059 buf_dp_size=buf_dp_size+size(paw_an_in1%k3xct1) 1060 end if 1061 if (paw_an_in1%has_vxcval==2) then 1062 buf_dp_size=buf_dp_size+size(paw_an_in1%vxc1_val) 1063 buf_dp_size=buf_dp_size+size(paw_an_in1%vxct1_val) 1064 end if 1065 if (paw_an_in1%has_vxc_ex==2) then 1066 buf_dp_size=buf_dp_size+size(paw_an_in1%vxc_ex) 1067 end if 1068 if (paw_an_in1%has_vhartree==2) then 1069 buf_dp_size=buf_dp_size+size(paw_an_in1%vh1) 1070 buf_dp_size=buf_dp_size+size(paw_an_in1%vht1) 1071 end if 1072 end do 1073 1074 !Fill in input buffers 1075 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 1076 LIBPAW_ALLOCATE(buf_dp ,(buf_dp_size)) 1077 indx_int=1;indx_dp=1 1078 do ij=1, my_natom 1079 paw_an_in1=>paw_an_in(ij) 1080 buf_int(indx_int)=my_atmtab(ij); indx_int=indx_int+1 1081 buf_int(indx_int)=paw_an_in1%itypat; indx_int=indx_int+1 1082 buf_int(indx_int)=paw_an_in1%nspden; indx_int=indx_int+1 1083 buf_int(indx_int)=paw_an_in1%cplex; indx_int=indx_int+1 1084 buf_int(indx_int)=paw_an_in1%mesh_size; indx_int=indx_int+1 1085 buf_int(indx_int)=paw_an_in1%angl_size; indx_int=indx_int+1 1086 buf_int(indx_int)=paw_an_in1%lm_size; indx_int=indx_int+1 1087 buf_int(indx_int)=paw_an_in1%nkxc1; indx_int=indx_int+1 1088 buf_int(indx_int)=paw_an_in1%nk3xc1; indx_int=indx_int+1 1089 buf_int(indx_int)=paw_an_in1%has_vxc; indx_int=indx_int+1 1090 buf_int(indx_int)=paw_an_in1%has_kxc; indx_int=indx_int+1 1091 buf_int(indx_int)=paw_an_in1%has_k3xc; indx_int=indx_int+1 1092 buf_int(indx_int)=paw_an_in1%has_vxcval; indx_int=indx_int+1 1093 buf_int(indx_int)=paw_an_in1%has_vxc_ex; indx_int=indx_int+1 1094 buf_int(indx_int)=paw_an_in1%has_vhartree; indx_int=indx_int+1 1095 v_size=0 1096 if (paw_an_in1%has_vxc>0) then 1097 v_size=size(paw_an_in1%vxc1,2) 1098 else if (paw_an_in1%has_kxc>0) then 1099 v_size=size(paw_an_in1%kxc1,2) 1100 else if (paw_an_in1%has_k3xc>0) then 1101 v_size=size(paw_an_in1%k3xc1,2) 1102 else if (paw_an_in1%has_vxcval>0) then 1103 v_size=size(paw_an_in1%vxc1_val,2) 1104 else if (paw_an_in1%has_vxc_ex>0) then 1105 v_size=size(paw_an_in1%vxc_ex,2) 1106 else if (paw_an_in1%has_vhartree>0) then 1107 v_size=size(paw_an_in1%vh1,2) 1108 end if 1109 buf_int(indx_int)=v_size;indx_int=indx_int+1 1110 if (allocated(paw_an_in1%lmselect)) then 1111 buf_int(indx_int)=1;indx_int=indx_int+1 1112 else 1113 buf_int(indx_int)=0;indx_int=indx_int+1 1114 end if 1115 nspden=paw_an_in1%nspden 1116 lm_size=paw_an_in1%lm_size 1117 cplx_mesh_size=paw_an_in1%cplex*paw_an_in1%mesh_size 1118 if (lm_size>0) then 1119 if (allocated(paw_an_in1%lmselect)) then 1120 do i1=1,lm_size 1121 if (paw_an_in1%lmselect(i1)) then 1122 buf_int(indx_int)=1 1123 else 1124 buf_int(indx_int)=0 1125 end if 1126 indx_int=indx_int+1 1127 end do 1128 end if 1129 end if 1130 if (paw_an_in1%has_vxc==2) then 1131 do i1=1,nspden 1132 do i2=1,v_size 1133 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vxc1(:,i2,i1) 1134 indx_dp=indx_dp+cplx_mesh_size 1135 end do 1136 end do 1137 do i1=1,nspden 1138 do i2=1,v_size 1139 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vxct1(:,i2,i1) 1140 indx_dp=indx_dp+cplx_mesh_size 1141 end do 1142 end do 1143 end if 1144 if (paw_an_in1%has_kxc==2.and.paw_an_in1%nkxc1>0) then 1145 do i1=1,paw_an_in1%nkxc1 1146 do i2=1,v_size 1147 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%kxc1(:,i2,i1) 1148 indx_dp=indx_dp+cplx_mesh_size 1149 end do 1150 end do 1151 do i1=1,paw_an_in1%nkxc1 1152 do i2=1,v_size 1153 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%kxct1(:,i2,i1) 1154 indx_dp=indx_dp+cplx_mesh_size 1155 end do 1156 end do 1157 end if 1158 if (paw_an_in1%has_k3xc==2.and.paw_an_in1%nk3xc1>0) then 1159 do i1=1,paw_an_in1%nk3xc1 1160 do i2=1,v_size 1161 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%k3xc1(:,i2,i1) 1162 indx_dp=indx_dp+cplx_mesh_size 1163 end do 1164 end do 1165 do i1=1,paw_an_in1%nk3xc1 1166 do i2=1,v_size 1167 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%k3xct1(:,i2,i1) 1168 indx_dp=indx_dp+cplx_mesh_size 1169 end do 1170 end do 1171 end if 1172 if (paw_an_in1%has_vxcval==2) then 1173 do i1=1,nspden 1174 do i2=1,v_size 1175 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vxc1_val(:,i2,i1) 1176 indx_dp=indx_dp+cplx_mesh_size 1177 end do 1178 end do 1179 do i1=1,nspden 1180 do i2=1,v_size 1181 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vxct1_val(:,i2,i1) 1182 indx_dp=indx_dp+cplx_mesh_size 1183 end do 1184 end do 1185 end if 1186 if (paw_an_in1%has_vxc_ex==2) then 1187 do i1=1,nspden 1188 do i2=1,v_size 1189 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vxc_ex(:,i2,i1) 1190 indx_dp=indx_dp+cplx_mesh_size 1191 end do 1192 end do 1193 end if 1194 if (paw_an_in1%has_vhartree==2) then 1195 do i1=1,nspden 1196 do i2=1,lm_size 1197 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vh1(:,i2,i1) 1198 indx_dp=indx_dp+cplx_mesh_size 1199 end do 1200 end do 1201 do i1=1,nspden 1202 do i2=1,lm_size 1203 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an_in1%vht1(:,i2,i1) 1204 indx_dp=indx_dp+cplx_mesh_size 1205 end do 1206 end do 1207 end if 1208 end do 1209 if (indx_int/=1+buf_int_size) then 1210 msg='Error (1) in paw_an_gather: wrong buffer sizes !' 1211 MSG_BUG(msg) 1212 end if 1213 if (indx_dp/=1+buf_dp_size) then 1214 msg='Error (2) in paw_an_gather: wrong buffer sizes !' 1215 MSG_BUG(msg) 1216 end if 1217 1218 !Communicate (1 gather for integers, 1 gather for reals) 1219 LIBPAW_ALLOCATE(count_int,(nproc_atom)) 1220 LIBPAW_ALLOCATE(displ_int,(nproc_atom)) 1221 LIBPAW_ALLOCATE(count_dp ,(nproc_atom)) 1222 LIBPAW_ALLOCATE(displ_dp ,(nproc_atom)) 1223 LIBPAW_ALLOCATE(count_tot,(2*nproc_atom)) 1224 bufsz(1)=buf_int_size; bufsz(2)=buf_dp_size 1225 call xmpi_allgather(bufsz,2,count_tot,comm_atom,ierr) 1226 do ij=1,nproc_atom 1227 count_int(ij)=count_tot(2*ij-1) 1228 count_dp (ij)=count_tot(2*ij) 1229 end do 1230 displ_int(1)=0;displ_dp(1)=0 1231 do ij=2,nproc_atom 1232 displ_int(ij)=displ_int(ij-1)+count_int(ij-1) 1233 displ_dp (ij)=displ_dp (ij-1)+count_dp (ij-1) 1234 end do 1235 buf_int_size_all=sum(count_int) 1236 buf_dp_size_all =sum(count_dp) 1237 LIBPAW_DEALLOCATE(count_tot) 1238 LIBPAW_ALLOCATE(buf_int_all,(buf_int_size_all)) 1239 LIBPAW_ALLOCATE(buf_dp_all ,(buf_dp_size_all)) 1240 call xmpi_allgatherv(buf_int,buf_int_size,buf_int_all,count_int,displ_int,comm_atom,ierr) 1241 call xmpi_allgatherv(buf_dp ,buf_dp_size ,buf_dp_all ,count_dp ,displ_dp ,comm_atom,ierr) 1242 LIBPAW_DEALLOCATE(count_int) 1243 LIBPAW_DEALLOCATE(displ_int) 1244 LIBPAW_DEALLOCATE(count_dp) 1245 LIBPAW_DEALLOCATE(displ_dp) 1246 1247 !Fill in output datastructure 1248 indx_int=1; indx_dp=1 1249 call paw_an_free(paw_an_gathered) 1250 call paw_an_nullify(paw_an_gathered) 1251 do iat=1,natom 1252 iatot=buf_int_all(indx_int); indx_int=indx_int+1 1253 paw_an_gathered1=>paw_an_gathered(iatot) 1254 paw_an_gathered1%itypat=buf_int_all(indx_int); indx_int=indx_int+1 1255 paw_an_gathered1%nspden=buf_int_all(indx_int); indx_int=indx_int+1 1256 paw_an_gathered1%cplex=buf_int_all(indx_int); indx_int=indx_int+1 1257 paw_an_gathered1%mesh_size=buf_int_all(indx_int); indx_int=indx_int+1 1258 paw_an_gathered1%angl_size=buf_int_all(indx_int); indx_int=indx_int+1 1259 paw_an_gathered1%lm_size=buf_int_all(indx_int); indx_int=indx_int+1 1260 paw_an_gathered1%nkxc1=buf_int_all(indx_int); indx_int=indx_int+1 1261 paw_an_gathered1%nk3xc1=buf_int_all(indx_int); indx_int=indx_int+1 1262 paw_an_gathered1%has_vxc=buf_int_all(indx_int); indx_int=indx_int+1 1263 paw_an_gathered1%has_kxc=buf_int_all(indx_int); indx_int=indx_int+1 1264 paw_an_gathered1%has_k3xc=buf_int_all(indx_int); indx_int=indx_int+1 1265 paw_an_gathered1%has_vxcval=buf_int_all(indx_int); indx_int=indx_int+1 1266 paw_an_gathered1%has_vxc_ex=buf_int_all(indx_int); indx_int=indx_int+1 1267 paw_an_gathered1%has_vhartree=buf_int_all(indx_int); indx_int=indx_int+1 1268 v_size=buf_int_all(indx_int); indx_int=indx_int+1 1269 has_lm_select=buf_int_all(indx_int); indx_int=indx_int+1 1270 nspden=paw_an_gathered1%nspden 1271 lm_size=paw_an_gathered1%lm_size 1272 nkxc1=paw_an_gathered1%nkxc1 1273 nk3xc1=paw_an_gathered1%nk3xc1 1274 cplx_mesh_size=paw_an_gathered1%cplex*paw_an_gathered1%mesh_size 1275 if (has_lm_select==1) then 1276 LIBPAW_ALLOCATE(paw_an_gathered1%lmselect,(lm_size)) 1277 if (lm_size>0) then 1278 do i1=1,lm_size 1279 if (buf_int_all(indx_int)==1) then 1280 paw_an_gathered1%lmselect(i1)=.TRUE.;indx_int=indx_int+1 1281 else 1282 paw_an_gathered1%lmselect(i1)=.FALSE.;indx_int=indx_int+1 1283 end if 1284 end do 1285 end if 1286 end if 1287 if (paw_an_gathered1%has_vxc>0) then 1288 LIBPAW_ALLOCATE(paw_an_gathered1%vxc1,(cplx_mesh_size,v_size,nspden)) 1289 LIBPAW_ALLOCATE(paw_an_gathered1%vxct1,(cplx_mesh_size,v_size,nspden)) 1290 if (paw_an_gathered1%has_vxc==2) then 1291 do i1=1,nspden 1292 do i2=1,v_size 1293 paw_an_gathered1%vxc1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1294 indx_dp=indx_dp+cplx_mesh_size 1295 end do 1296 end do 1297 do i1=1,nspden 1298 do i2=1,v_size 1299 paw_an_gathered1%vxct1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1300 indx_dp=indx_dp+cplx_mesh_size 1301 end do 1302 end do 1303 end if 1304 end if 1305 if (paw_an_gathered1%has_kxc>0) then 1306 LIBPAW_ALLOCATE(paw_an_gathered1%kxc1,(cplx_mesh_size,v_size,nkxc1)) 1307 LIBPAW_ALLOCATE(paw_an_gathered1%kxct1,(cplx_mesh_size,v_size,nkxc1)) 1308 if (paw_an_gathered1%has_kxc==2.and.nkxc1>0) then 1309 do i1=1,nkxc1 1310 do i2=1,v_size 1311 paw_an_gathered1%kxc1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1312 indx_dp=indx_dp+cplx_mesh_size 1313 end do 1314 end do 1315 do i1=1,nkxc1 1316 do i2=1,v_size 1317 paw_an_gathered1%kxct1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1318 indx_dp=indx_dp+cplx_mesh_size 1319 end do 1320 end do 1321 end if 1322 end if 1323 if (paw_an_gathered1%has_k3xc>0) then 1324 LIBPAW_ALLOCATE(paw_an_gathered1%k3xc1,(cplx_mesh_size,v_size,nk3xc1)) 1325 LIBPAW_ALLOCATE(paw_an_gathered1%k3xct1,(cplx_mesh_size,v_size,nk3xc1)) 1326 if (paw_an_gathered1%has_k3xc==2.and.nk3xc1>0) then 1327 do i1=1,nk3xc1 1328 do i2=1,v_size 1329 paw_an_gathered1%k3xc1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1330 indx_dp=indx_dp+cplx_mesh_size 1331 end do 1332 end do 1333 do i1=1,nk3xc1 1334 do i2=1,v_size 1335 paw_an_gathered1%k3xct1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1336 indx_dp=indx_dp+cplx_mesh_size 1337 end do 1338 end do 1339 end if 1340 end if 1341 if (paw_an_gathered1%has_vxcval>0) then 1342 LIBPAW_ALLOCATE(paw_an_gathered1%vxc1_val,(cplx_mesh_size,v_size,nspden)) 1343 LIBPAW_ALLOCATE(paw_an_gathered1%vxct1_val,(cplx_mesh_size,v_size,nspden)) 1344 if (paw_an_gathered1%has_vxcval==2) then 1345 do i1=1,nspden 1346 do i2=1,v_size 1347 paw_an_gathered1%vxc1_val(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1348 indx_dp=indx_dp+cplx_mesh_size 1349 end do 1350 end do 1351 do i1=1,nspden 1352 do i2=1,v_size 1353 paw_an_gathered1%vxct1_val(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1354 indx_dp=indx_dp+cplx_mesh_size 1355 end do 1356 end do 1357 end if 1358 end if 1359 if (paw_an_gathered1%has_vxc_ex>0) then 1360 LIBPAW_ALLOCATE(paw_an_gathered1%vxc_ex,(cplx_mesh_size,v_size,nspden)) 1361 if (paw_an_gathered1%has_vxc_ex==2) then 1362 do i1=1,nspden 1363 do i2=1,v_size 1364 paw_an_gathered1%vxc_ex(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1365 indx_dp=indx_dp+cplx_mesh_size 1366 end do 1367 end do 1368 end if 1369 end if 1370 if (paw_an_gathered1%has_vhartree>0) then 1371 LIBPAW_ALLOCATE(paw_an_gathered1%vh1,(cplx_mesh_size,lm_size,nspden)) 1372 LIBPAW_ALLOCATE(paw_an_gathered1%vht1,(cplx_mesh_size,lm_size,nspden)) 1373 if (paw_an_gathered1%has_vhartree==2) then 1374 do i1=1,nspden 1375 do i2=1,lm_size 1376 paw_an_gathered1%vh1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1377 indx_dp=indx_dp+cplx_mesh_size 1378 end do 1379 end do 1380 do i1=1,nspden 1381 do i2=1,lm_size 1382 paw_an_gathered1%vht1(:,i2,i1)=buf_dp_all(indx_dp:indx_dp+cplx_mesh_size-1) 1383 indx_dp=indx_dp+cplx_mesh_size 1384 end do 1385 end do 1386 end if 1387 end if 1388 end do ! iat 1389 1390 !Free buffers 1391 LIBPAW_DEALLOCATE(buf_int) 1392 LIBPAW_DEALLOCATE(buf_int_all) 1393 LIBPAW_DEALLOCATE(buf_dp) 1394 LIBPAW_DEALLOCATE(buf_dp_all) 1395 1396 !Destroy atom table 1397 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 1398 1399 end subroutine paw_an_gather
m_paw_an/paw_an_init [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_init
FUNCTION
Initialize a paw_an data type.
INPUTS
mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc comm_atom=--optional-- MPI communicator over atoms
SIDE EFFECTS
Paw_an(:)<type(paw_an_type)>=PAW arrays given on ANgular mesh or ANgular moments. Initialized in output
PARENTS
bethe_salpeter,dfpt_scfcv,paw_qpscgw,respfn,scfcv,screening,sigma
CHILDREN
SOURCE
228 subroutine paw_an_init(Paw_an,natom,ntypat,nkxc1,nk3xc1,nspden,cplex,pawxcdev,typat,Pawang,Pawtab,& 229 & has_vhartree,has_vxc,has_vxcval,has_kxc,has_k3xc,has_vxc_ex, & ! optional arguments 230 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 231 232 233 !This section has been created automatically by the script Abilint (TD). 234 !Do not modify the following lines by hand. 235 #undef ABI_FUNC 236 #define ABI_FUNC 'paw_an_init' 237 !End of the abilint section 238 239 implicit none 240 241 !Arguments ------------------------------------ 242 !scalars 243 integer,intent(in) :: natom,nkxc1,nk3xc1,ntypat,cplex,nspden,pawxcdev 244 integer,optional,intent(in) :: has_vhartree,has_vxc,has_vxcval,has_kxc,has_k3xc,has_vxc_ex 245 integer,optional,intent(in) :: comm_atom 246 !arrays 247 integer,intent(in) :: typat(natom) 248 integer,optional,target,intent(in) :: mpi_atmtab(:) 249 type(Pawang_type),intent(in) :: Pawang 250 type(Pawtab_type),intent(in) :: Pawtab(ntypat) 251 type(Paw_an_type),intent(inout) :: Paw_an(:) 252 253 !Local variables------------------------------- 254 !scalars 255 integer :: iat,iat1,itypat,lm_size,my_comm_atom,my_natom,v_size 256 logical :: my_atmtab_allocated,paral_atom 257 !arrays 258 integer,pointer :: my_atmtab(:) 259 260 ! ************************************************************************* 261 262 !@Paw_an_type 263 264 !Set up parallelism over atoms 265 my_natom=size(Paw_an);if (my_natom==0) return 266 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 267 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 268 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 269 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 270 271 do iat=1,my_natom 272 iat1=iat;if (paral_atom) iat1=my_atmtab(iat) 273 itypat=typat(iat1) 274 275 lm_size =Pawtab(itypat)%lcut_size**2 276 Paw_an(iat)%angl_size =Pawang%angl_size 277 Paw_an(iat)%cplex =cplex 278 Paw_an(iat)%itypat =itypat 279 Paw_an(iat)%lm_size =lm_size 280 Paw_an(iat)%mesh_size =Pawtab(itypat)%mesh_size 281 Paw_an(iat)%nkxc1 =nkxc1 282 Paw_an(iat)%nk3xc1 =nk3xc1 283 Paw_an(iat)%nspden =nspden 284 285 ! === Non-zero LM-moments of "one-center" densities/potentials === 286 ! * Filled in pawdenpot. 287 LIBPAW_ALLOCATE(Paw_an(iat)%lmselect,(lm_size)) 288 289 v_size=Paw_an(iat)%lm_size ; if (pawxcdev==0) v_size=Paw_an(iat)%angl_size 290 291 ! === XC potential inside the sphere === 292 ! * LM-moments of potential if pawxcdev/=0 293 ! * (theta,phi) values of potential if pawxcdev=0 294 Paw_an(iat)%has_vxc=0 295 if (PRESENT(has_vxc)) then 296 if (has_vxc>0) then 297 Paw_an(iat)%has_vxc=1 298 LIBPAW_ALLOCATE(Paw_an(iat)%vxc1 ,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 299 LIBPAW_ALLOCATE(Paw_an(iat)%vxct1,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 300 Paw_an(iat)%vxc1=zero;Paw_an(iat)%vxct1=zero 301 end if 302 end if 303 304 ! ========================== 305 ! === Optional arguments === 306 ! ========================== 307 308 ! * XC potential inside PAW spheres generated by valence electrons 309 Paw_an(iat)%has_vxcval=0 310 if (PRESENT(has_vxcval)) then 311 if (has_vxcval>0) then 312 Paw_an(iat)%has_vxcval=1 313 LIBPAW_ALLOCATE(Paw_an(iat)%vxc1_val ,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 314 LIBPAW_ALLOCATE(Paw_an(iat)%vxct1_val,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 315 Paw_an(iat)%vxc1_val=zero;Paw_an(iat)%vxct1_val=zero 316 end if 317 end if 318 319 ! * Hartree potential LM-moments inside the sphere 320 Paw_an(iat)%has_vhartree=0 321 if (PRESENT(has_vhartree)) then 322 if (has_vhartree>0) then 323 Paw_an(iat)%has_vhartree=1 324 LIBPAW_ALLOCATE(Paw_an(iat)%vh1,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 325 LIBPAW_ALLOCATE(Paw_an(iat)%vht1,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 326 Paw_an(iat)%vh1=zero;Paw_an(iat)%vht1=zero 327 end if 328 end if 329 330 ! xc kernels inside the sphere 331 Paw_an(iat)%has_kxc=0 332 if (PRESENT(has_kxc)) then 333 if (has_kxc>0) then 334 Paw_an(iat)%has_kxc=1 335 LIBPAW_ALLOCATE(Paw_an(iat)%kxc1 ,(cplex*Paw_an(iat)%mesh_size,v_size,nkxc1)) 336 LIBPAW_ALLOCATE(Paw_an(iat)%kxct1,(cplex*Paw_an(iat)%mesh_size,v_size,nkxc1)) 337 if (nkxc1>0) then 338 Paw_an(iat)%kxc1=zero;Paw_an(iat)%kxct1=zero 339 end if 340 end if 341 end if 342 343 ! xc kernel derivatives inside the sphere 344 Paw_an(iat)%has_k3xc=0 345 if (PRESENT(has_k3xc)) then 346 if (has_k3xc>0) then 347 Paw_an(iat)%has_k3xc=1 348 LIBPAW_ALLOCATE(Paw_an(iat)%k3xc1 ,(cplex*Paw_an(iat)%mesh_size,v_size,nk3xc1)) 349 LIBPAW_ALLOCATE(Paw_an(iat)%k3xct1,(cplex*Paw_an(iat)%mesh_size,v_size,nk3xc1)) 350 if (nk3xc1>0) then 351 Paw_an(iat)%k3xc1=zero;Paw_an(iat)%k3xct1=zero 352 end if 353 end if 354 end if 355 356 ! local exact-exchange potential inside the sphere 357 Paw_an(iat)%has_vxc_ex=0 358 if (PRESENT(has_vxc_ex)) then 359 if (has_vxc_ex>0.and.Pawtab(itypat)%useexexch>0) then 360 Paw_an(iat)%has_vxc_ex=1 361 LIBPAW_ALLOCATE(Paw_an(iat)%vxc_ex,(cplex*Paw_an(iat)%mesh_size,v_size,nspden)) 362 Paw_an(iat)%vxc_ex=zero 363 end if 364 end if 365 366 end do !iat 367 368 !Destroy atom table used for parallelism 369 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 370 371 end subroutine paw_an_init
m_paw_an/paw_an_isendreceive_fillbuffer [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_isendreceive_fillbuffer
FUNCTION
Extract from paw_an 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_an_send= number of sent atoms paw_an= data structure from which are extract buffer int and buffer dp
OUTPUT
buf_int : buffer of integers to be send in a send operation buf_int_size : size of buffer of integers to be send in a send operation buf_dp : buffer of double precision numbers to be send in a send operation buf_dp_size : size of buffer of double precision numbers to be send in a send operation
PARENTS
m_paw_an
CHILDREN
SOURCE
2061 subroutine paw_an_isendreceive_fillbuffer(paw_an, atmtab_send,atm_indx_send,npaw_an_send,& 2062 & buf_int,buf_int_size,buf_dp,buf_dp_size) 2063 2064 2065 !This section has been created automatically by the script Abilint (TD). 2066 !Do not modify the following lines by hand. 2067 #undef ABI_FUNC 2068 #define ABI_FUNC 'paw_an_isendreceive_fillbuffer' 2069 !End of the abilint section 2070 2071 implicit none 2072 !Arguments ------------------------------------ 2073 !scalars 2074 integer,intent(out) :: buf_int_size,buf_dp_size 2075 integer,intent(in) :: npaw_an_send 2076 !arrays 2077 integer,intent(in) :: atmtab_send(:),atm_indx_send(:) 2078 integer,allocatable,intent(out) :: buf_int(:) 2079 real(dp),allocatable,intent(out):: buf_dp(:) 2080 type(paw_an_type),target,intent(in) :: paw_an(:) 2081 2082 !Local variables------------------------------- 2083 !scalars 2084 integer :: cplx_mesh_size,i1,i2,iatom_tot,ij,indx_int,indx_dp 2085 integer :: ipaw_an_send,lm_size,nspden,v_size 2086 character(len=500) :: msg 2087 type(Paw_an_type),pointer :: paw_an1 2088 !arrays 2089 2090 ! ********************************************************************* 2091 2092 !Compute sizes of buffers 2093 buf_int_size=0 ; buf_dp_size=0 2094 do ipaw_an_send=1,npaw_an_send 2095 iatom_tot=atmtab_send(ipaw_an_send) 2096 ij = atm_indx_send(iatom_tot) 2097 paw_an1=>paw_an(ij) 2098 buf_int_size=buf_int_size+17+size(paw_an1%lmselect) 2099 if (paw_an1%has_vxc==2) then 2100 buf_dp_size=buf_dp_size+size(paw_an1%vxc1) 2101 buf_dp_size=buf_dp_size+size(paw_an1%vxct1) 2102 end if 2103 if (paw_an1%has_kxc==2) then 2104 buf_dp_size=buf_dp_size+size(paw_an1%kxc1) 2105 buf_dp_size=buf_dp_size+size(paw_an1%kxct1) 2106 end if 2107 if (paw_an1%has_k3xc==2) then 2108 buf_dp_size=buf_dp_size+size(paw_an1%k3xc1) 2109 buf_dp_size=buf_dp_size+size(paw_an1%k3xct1) 2110 end if 2111 if (paw_an1%has_vxcval==2) then 2112 buf_dp_size=buf_dp_size+size(paw_an1%vxc1_val) 2113 buf_dp_size=buf_dp_size+size(paw_an1%vxct1_val) 2114 end if 2115 if (paw_an1%has_vxc_ex==2) then 2116 buf_dp_size=buf_dp_size+size(paw_an1%vxc_ex) 2117 end if 2118 if (paw_an1%has_vhartree==2) then 2119 buf_dp_size=buf_dp_size+size(paw_an1%vh1) 2120 buf_dp_size=buf_dp_size+size(paw_an1%vht1) 2121 end if 2122 end do 2123 2124 !Fill in input buffers 2125 LIBPAW_ALLOCATE(buf_int,(buf_int_size)) 2126 LIBPAW_ALLOCATE(buf_dp ,(buf_dp_size)) 2127 indx_int=1;indx_dp=1 2128 do ipaw_an_send=1,npaw_an_send 2129 iatom_tot=atmtab_send(ipaw_an_send) 2130 ij = atm_indx_send(iatom_tot) 2131 paw_an1=>paw_an(ij) 2132 buf_int(indx_int)=iatom_tot; indx_int=indx_int+1 2133 buf_int(indx_int)=paw_an1%itypat; indx_int=indx_int+1 2134 buf_int(indx_int)=paw_an1%nspden; indx_int=indx_int+1 2135 buf_int(indx_int)=paw_an1%cplex; indx_int=indx_int+1 2136 buf_int(indx_int)=paw_an1%mesh_size; indx_int=indx_int+1 2137 buf_int(indx_int)=paw_an1%angl_size; indx_int=indx_int+1 2138 buf_int(indx_int)=paw_an1%lm_size; indx_int=indx_int+1 2139 buf_int(indx_int)=paw_an1%nkxc1; indx_int=indx_int+1 2140 buf_int(indx_int)=paw_an1%nk3xc1; indx_int=indx_int+1 2141 buf_int(indx_int)=paw_an1%has_vxc; indx_int=indx_int+1 2142 buf_int(indx_int)=paw_an1%has_kxc; indx_int=indx_int+1 2143 buf_int(indx_int)=paw_an1%has_k3xc; indx_int=indx_int+1 2144 buf_int(indx_int)=paw_an1%has_vxcval; indx_int=indx_int+1 2145 buf_int(indx_int)=paw_an1%has_vxc_ex; indx_int=indx_int+1 2146 buf_int(indx_int)=paw_an1%has_vhartree; indx_int=indx_int+1 2147 v_size=0 2148 if (paw_an1%has_vxc>0) then 2149 v_size=size(paw_an1%vxc1,2) 2150 else if (paw_an1%has_kxc>0) then 2151 v_size=size(paw_an1%kxc1,2) 2152 else if (paw_an1%has_k3xc>0) then 2153 v_size=size(paw_an1%k3xc1,2) 2154 else if (paw_an1%has_vxcval>0) then 2155 v_size=size(paw_an1%vxc1_val,2) 2156 else if (paw_an1%has_vxc_ex>0) then 2157 v_size=size(paw_an1%vxc_ex,2) 2158 else if (paw_an1%has_vhartree>0) then 2159 v_size=size(paw_an1%vh1,2) 2160 end if 2161 buf_int(indx_int)=v_size;indx_int=indx_int+1 2162 if (allocated(paw_an1%lmselect)) then 2163 buf_int(indx_int)=1;indx_int=indx_int+1 2164 else 2165 buf_int(indx_int)=0;indx_int=indx_int+1 2166 end if 2167 nspden=paw_an1%nspden 2168 lm_size=paw_an1%lm_size 2169 cplx_mesh_size=paw_an1%cplex*paw_an1%mesh_size 2170 if (lm_size>0) then 2171 if (allocated(paw_an1%lmselect)) then 2172 do i1=1,lm_size 2173 if (paw_an1%lmselect(i1)) then 2174 buf_int(indx_int)=1 2175 else 2176 buf_int(indx_int)=0 2177 end if 2178 indx_int=indx_int+1 2179 end do 2180 end if 2181 end if 2182 if (paw_an1%has_vxc==2) then 2183 do i1=1,nspden 2184 do i2=1,v_size 2185 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vxc1(:,i2,i1) 2186 indx_dp=indx_dp+cplx_mesh_size 2187 end do 2188 end do 2189 do i1=1,nspden 2190 do i2=1,v_size 2191 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vxct1(:,i2,i1) 2192 indx_dp=indx_dp+cplx_mesh_size 2193 end do 2194 end do 2195 end if 2196 if (paw_an1%has_kxc==2.and.paw_an1%nkxc1>0) then 2197 do i1=1,paw_an1%nkxc1 2198 do i2=1,v_size 2199 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%kxc1(:,i2,i1) 2200 indx_dp=indx_dp+cplx_mesh_size 2201 end do 2202 end do 2203 do i1=1,paw_an1%nkxc1 2204 do i2=1,v_size 2205 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%kxct1(:,i2,i1) 2206 indx_dp=indx_dp+cplx_mesh_size 2207 end do 2208 end do 2209 end if 2210 if (paw_an1%has_k3xc==2.and.paw_an1%nk3xc1>0) then 2211 do i1=1,paw_an1%nk3xc1 2212 do i2=1,v_size 2213 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%k3xc1(:,i2,i1) 2214 indx_dp=indx_dp+cplx_mesh_size 2215 end do 2216 end do 2217 do i1=1,paw_an1%nk3xc1 2218 do i2=1,v_size 2219 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%k3xct1(:,i2,i1) 2220 indx_dp=indx_dp+cplx_mesh_size 2221 end do 2222 end do 2223 end if 2224 2225 if (paw_an1%has_vxcval==2) then 2226 do i1=1,nspden 2227 do i2=1,v_size 2228 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vxc1_val(:,i2,i1) 2229 indx_dp=indx_dp+cplx_mesh_size 2230 end do 2231 end do 2232 do i1=1,nspden 2233 do i2=1,v_size 2234 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vxct1_val(:,i2,i1) 2235 indx_dp=indx_dp+cplx_mesh_size 2236 end do 2237 end do 2238 end if 2239 if (paw_an1%has_vxc_ex==2) then 2240 do i1=1,nspden 2241 do i2=1,v_size 2242 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vxc_ex(:,i2,i1) 2243 indx_dp=indx_dp+cplx_mesh_size 2244 end do 2245 end do 2246 end if 2247 if (paw_an1%has_vhartree==2) then 2248 do i1=1,nspden 2249 do i2=1,lm_size 2250 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vh1(:,i2,i1) 2251 indx_dp=indx_dp+cplx_mesh_size 2252 end do 2253 end do 2254 do i1=1,nspden 2255 do i2=1,lm_size 2256 buf_dp(indx_dp:indx_dp+cplx_mesh_size-1)=paw_an1%vht1(:,i2,i1) 2257 indx_dp=indx_dp+cplx_mesh_size 2258 end do 2259 end do 2260 end if 2261 end do 2262 if ((indx_int-1/=buf_int_size).or.(indx_dp-1/=buf_dp_size)) then 2263 write(msg,'(4(a,i10))') 'Wrong buffer sizes: buf_int =',buf_int_size,'/',indx_int-1,& 2264 & ' buf_dp =',buf_dp_size ,'/',indx_dp-1 2265 MSG_BUG(msg) 2266 end if 2267 2268 end subroutine paw_an_isendreceive_fillbuffer
m_paw_an/paw_an_isendreceive_getbuffer [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_isendreceive_getbuffer
FUNCTION
Fill a paw_an structure with the buffers received in a receive operation This buffer should have been first extracted by a call to paw_an_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_an_send= number of sent atoms
OUTPUT
paw_an= output datastructure filled with buffers receive in a receive operation
PARENTS
m_paw_an
CHILDREN
SOURCE
1849 subroutine paw_an_isendreceive_getbuffer(paw_an,npaw_an_send,atm_indx_recv,buf_int,buf_dp) 1850 1851 1852 !This section has been created automatically by the script Abilint (TD). 1853 !Do not modify the following lines by hand. 1854 #undef ABI_FUNC 1855 #define ABI_FUNC 'paw_an_isendreceive_getbuffer' 1856 !End of the abilint section 1857 1858 implicit none 1859 1860 !Arguments ------------------------------------ 1861 !scalars 1862 integer,intent(in) :: npaw_an_send 1863 !arrays 1864 integer,intent(in) ::atm_indx_recv(:),buf_int(:) 1865 real(dp),intent(in):: buf_dp(:) 1866 type(paw_an_type),target,intent(inout) :: paw_an(:) 1867 1868 !Local variables------------------------------- 1869 !scalars 1870 integer :: buf_int_size,buf_dp_size,cplx_mesh_size,has_lm_select,i1,i2 1871 integer :: iat,iatot,ij,indx_int,indx_dp,lm_size,nkxc1,nk3xc1,nspden,v_size 1872 character(len=500) :: msg 1873 type(paw_an_type),pointer :: paw_an1 1874 !arrays 1875 1876 ! ********************************************************************* 1877 1878 buf_int_size=size(buf_int) 1879 buf_dp_size=size(buf_dp) 1880 indx_int=1; indx_dp=1 1881 1882 do ij=1,npaw_an_send 1883 iatot=buf_int(indx_int); indx_int=indx_int+1 1884 iat= atm_indx_recv(iatot) 1885 paw_an1=>paw_an(iat) 1886 paw_an1%itypat=buf_int(indx_int); indx_int=indx_int+1 1887 paw_an1%nspden=buf_int(indx_int); indx_int=indx_int+1 1888 paw_an1%cplex=buf_int(indx_int); indx_int=indx_int+1 1889 paw_an1%mesh_size=buf_int(indx_int); indx_int=indx_int+1 1890 paw_an1%angl_size=buf_int(indx_int); indx_int=indx_int+1 1891 paw_an1%lm_size=buf_int(indx_int); indx_int=indx_int+1 1892 paw_an1%nkxc1=buf_int(indx_int); indx_int=indx_int+1 1893 paw_an1%nk3xc1=buf_int(indx_int); indx_int=indx_int+1 1894 paw_an1%has_vxc=buf_int(indx_int); indx_int=indx_int+1 1895 paw_an1%has_kxc=buf_int(indx_int); indx_int=indx_int+1 1896 paw_an1%has_k3xc=buf_int(indx_int); indx_int=indx_int+1 1897 paw_an1%has_vxcval=buf_int(indx_int); indx_int=indx_int+1 1898 paw_an1%has_vxc_ex=buf_int(indx_int); indx_int=indx_int+1 1899 paw_an1%has_vhartree=buf_int(indx_int); indx_int=indx_int+1 1900 v_size=buf_int(indx_int); indx_int=indx_int+1 1901 has_lm_select=buf_int(indx_int); indx_int=indx_int+1 1902 nspden=paw_an1%nspden 1903 lm_size=paw_an1%lm_size 1904 nkxc1=paw_an1%nkxc1 1905 nk3xc1=paw_an1%nk3xc1 1906 cplx_mesh_size=paw_an1%cplex*paw_an1%mesh_size 1907 if (has_lm_select==1) then 1908 LIBPAW_ALLOCATE(paw_an1%lmselect,(lm_size)) 1909 if (lm_size>0) then 1910 do i1=1,lm_size 1911 if (buf_int(indx_int)==1) then 1912 paw_an1%lmselect(i1)=.TRUE.;indx_int=indx_int+1 1913 else 1914 paw_an1%lmselect(i1)=.FALSE.;indx_int=indx_int+1 1915 end if 1916 end do 1917 end if 1918 end if 1919 if (paw_an1%has_vxc>0) then 1920 LIBPAW_ALLOCATE(paw_an1%vxc1,(cplx_mesh_size,v_size,nspden)) 1921 LIBPAW_ALLOCATE(paw_an1%vxct1,(cplx_mesh_size,v_size,nspden)) 1922 if (paw_an1%has_vxc==2) then 1923 do i1=1,nspden 1924 do i2=1,v_size 1925 paw_an1%vxc1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1926 indx_dp=indx_dp+cplx_mesh_size 1927 end do 1928 end do 1929 do i1=1,nspden 1930 do i2=1,v_size 1931 paw_an1%vxct1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1932 indx_dp=indx_dp+cplx_mesh_size 1933 end do 1934 end do 1935 end if 1936 end if 1937 if (paw_an1%has_kxc>0) then 1938 LIBPAW_ALLOCATE(paw_an1%kxc1,(cplx_mesh_size,v_size,nkxc1)) 1939 LIBPAW_ALLOCATE(paw_an1%kxct1,(cplx_mesh_size,v_size,nkxc1)) 1940 if (paw_an1%has_kxc==2.and.nkxc1>0) then 1941 do i1=1,nkxc1 1942 do i2=1,v_size 1943 paw_an1%kxc1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1944 indx_dp=indx_dp+cplx_mesh_size 1945 end do 1946 end do 1947 do i1=1,nkxc1 1948 do i2=1,v_size 1949 paw_an1%kxct1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1950 indx_dp=indx_dp+cplx_mesh_size 1951 end do 1952 end do 1953 end if 1954 end if 1955 if (paw_an1%has_k3xc>0) then 1956 LIBPAW_ALLOCATE(paw_an1%k3xc1,(cplx_mesh_size,v_size,nk3xc1)) 1957 LIBPAW_ALLOCATE(paw_an1%k3xct1,(cplx_mesh_size,v_size,nk3xc1)) 1958 if (paw_an1%has_k3xc==2.and.nk3xc1>0) then 1959 do i1=1,nk3xc1 1960 do i2=1,v_size 1961 paw_an1%k3xc1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1962 indx_dp=indx_dp+cplx_mesh_size 1963 end do 1964 end do 1965 do i1=1,nk3xc1 1966 do i2=1,v_size 1967 paw_an1%k3xct1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1968 indx_dp=indx_dp+cplx_mesh_size 1969 end do 1970 end do 1971 end if 1972 end if 1973 if (paw_an1%has_vxcval>0) then 1974 LIBPAW_ALLOCATE(paw_an1%vxc1_val,(cplx_mesh_size,v_size,nspden)) 1975 LIBPAW_ALLOCATE(paw_an1%vxct1_val,(cplx_mesh_size,v_size,nspden)) 1976 if (paw_an1%has_vxcval==2) then 1977 do i1=1,nspden 1978 do i2=1,v_size 1979 paw_an1%vxc1_val(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1980 indx_dp=indx_dp+cplx_mesh_size 1981 end do 1982 end do 1983 do i1=1,nspden 1984 do i2=1,v_size 1985 paw_an1%vxct1_val(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1986 indx_dp=indx_dp+cplx_mesh_size 1987 end do 1988 end do 1989 end if 1990 end if 1991 if (paw_an1%has_vxc_ex>0) then 1992 LIBPAW_ALLOCATE(paw_an1%vxc_ex,(cplx_mesh_size,v_size,nspden)) 1993 if (paw_an1%has_vxc_ex==2) then 1994 do i1=1,nspden 1995 do i2=1,v_size 1996 paw_an1%vxc_ex(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 1997 indx_dp=indx_dp+cplx_mesh_size 1998 end do 1999 end do 2000 end if 2001 end if 2002 if (paw_an1%has_vhartree>0) then 2003 LIBPAW_ALLOCATE(paw_an1%vh1,(cplx_mesh_size,lm_size,nspden)) 2004 LIBPAW_ALLOCATE(paw_an1%vht1,(cplx_mesh_size,lm_size,nspden)) 2005 if (paw_an1%has_vhartree==2) then 2006 do i1=1,nspden 2007 do i2=1,lm_size 2008 paw_an1%vh1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 2009 indx_dp=indx_dp+cplx_mesh_size 2010 end do 2011 end do 2012 do i1=1,nspden 2013 do i2=1,lm_size 2014 paw_an1%vht1(:,i2,i1)=buf_dp(indx_dp:indx_dp+cplx_mesh_size-1) 2015 indx_dp=indx_dp+cplx_mesh_size 2016 end do 2017 end do 2018 end if 2019 end if 2020 end do ! iat 2021 if ((indx_int/=1+buf_int_size).or.(indx_dp/=1+buf_dp_size)) then 2022 write(msg,'(a,i10,a,i10)') 'Wrong buffer sizes: buf_int_size=',buf_int_size,' buf_dp_size=',buf_dp_size 2023 MSG_BUG(msg) 2024 end if 2025 2026 end subroutine paw_an_isendreceive_getbuffer
m_paw_an/paw_an_nullify [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_nullify
FUNCTION
Nullify pointers and flags in a paw_an structure
INPUTS
SIDE EFFECTS
Paw_an(:)<type(paw_an_type)>=PAW arrays given on ANgular mesh or ANgular moments. Nullified in output
PARENTS
bethe_salpeter,dfpt_scfcv,m_paw_an,paw_qpscgw,respfn,scfcv,screening sigma
CHILDREN
SOURCE
493 subroutine paw_an_nullify(Paw_an) 494 495 496 !This section has been created automatically by the script Abilint (TD). 497 !Do not modify the following lines by hand. 498 #undef ABI_FUNC 499 #define ABI_FUNC 'paw_an_nullify' 500 !End of the abilint section 501 502 implicit none 503 504 !Arguments ------------------------------------ 505 type(Paw_an_type),intent(inout) :: Paw_an(:) 506 507 !Local variables------------------------------- 508 integer :: iat,natom 509 510 ! ************************************************************************* 511 512 !@Paw_an_type 513 ! MGPAW: This one could be removed/renamed, 514 ! variables can be initialized in the datatype declaration 515 ! Do we need to expose this in the public API? 516 517 natom=SIZE(Paw_an(:));if (natom==0) return 518 519 do iat=1,natom 520 ! Set all has_* flags to zero. 521 Paw_an(iat)%has_kxc =0 522 Paw_an(iat)%has_k3xc =0 523 Paw_an(iat)%has_vhartree =0 524 Paw_an(iat)%has_vxc =0 525 Paw_an(iat)%has_vxcval =0 526 Paw_an(iat)%has_vxc_ex =0 527 end do 528 529 end subroutine paw_an_nullify
m_paw_an/paw_an_print [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_print
FUNCTION
Reports basic info on a paw_an datastructure
INPUTS
[unit]=the unit number for output [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_an is distributed, natom is different from size(Paw_an).
OUTPUT
(only writing)
PARENTS
CHILDREN
SOURCE
766 subroutine paw_an_print(Paw_an,unit,mode_paral, & 767 & mpi_atmtab,comm_atom,natom) 768 769 770 !This section has been created automatically by the script Abilint (TD). 771 !Do not modify the following lines by hand. 772 #undef ABI_FUNC 773 #define ABI_FUNC 'paw_an_print' 774 !End of the abilint section 775 776 implicit none 777 778 !Arguments ------------------------------------ 779 !scalars 780 integer,optional,intent(in) :: comm_atom,natom,unit 781 character(len=4),optional,intent(in) :: mode_paral 782 !arrays 783 integer,optional,target,intent(in) :: mpi_atmtab(:) 784 type(Paw_an_type),intent(in) :: Paw_an(:) 785 786 !Local variables------------------------------- 787 !scalars 788 integer :: iatom,iatom_tot,my_comm_atom,my_natom,my_unt,size_paw_an 789 logical :: my_atmtab_allocated,paral_atom 790 character(len=4) :: my_mode 791 character(len=500) :: msg 792 !arrays 793 integer,pointer :: my_atmtab(:) 794 795 ! ************************************************************************* 796 797 !@Paw_an_type 798 799 size_paw_an=SIZE(Paw_an) 800 my_unt =std_out; if (PRESENT(unit )) my_unt =unit 801 my_mode ='PERS' ; if (PRESENT(mode_paral)) my_mode =mode_paral 802 my_natom=size_paw_an; if (PRESENT(natom)) my_natom=natom 803 804 !Set up parallelism over atoms 805 paral_atom=(present(comm_atom).and.my_natom/=size_paw_an) 806 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 807 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 808 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,my_natom,my_natom_ref=size_paw_an) 809 810 write(msg,'(3a)')ch10,' === Content of the pawfgrtab datatype === ',ch10 811 call wrtout(my_unt,msg,my_mode) 812 813 do iatom=1,my_natom 814 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 815 write(msg,'(a)')' ' 816 call wrtout(my_unt,msg,my_mode) 817 write(msg,'(a,i4)')' ****************************** iatom= ' , iatom_tot 818 call wrtout(my_unt,msg,my_mode) 819 write(msg,'(a,i4)')' Dimension of paw angular mesh= ',paw_an(iatom)%angl_size 820 call wrtout(my_unt,msg,my_mode) 821 write(msg,'(a,i4)')' cplex (1 if potentials/densities are real, 2 if they are complex)= ',& 822 & paw_an(iatom)%cplex 823 call wrtout(my_unt,msg,my_mode) 824 write(msg,'(a,i4)')' has_kxc = ',paw_an(iatom)%has_kxc 825 call wrtout(my_unt,msg,my_mode) 826 write(msg,'(a,i4)')' has_k3xc = ',paw_an(iatom)%has_k3xc 827 call wrtout(my_unt,msg,my_mode) 828 write(msg,'(a,i4)')' has_vhartree= ',paw_an(iatom)%has_vhartree 829 call wrtout(my_unt,msg,my_mode) 830 write(msg,'(a,i4)')' has_vxc = ',paw_an(iatom)%has_vxc 831 call wrtout(my_unt,msg,my_mode) 832 write(msg,'(a,i4)')' has_vxcval = ',paw_an(iatom)%has_vxcval 833 call wrtout(my_unt,msg,my_mode) 834 write(msg,'(a,i4)')' has_vxc_ex = ',paw_an(iatom)%has_vxc_ex 835 call wrtout(my_unt,msg,my_mode) 836 write(msg,'(a,i4)')' Atome type = ',paw_an(iatom)%itypat 837 call wrtout(my_unt,msg,my_mode) 838 write(msg,'(a,i4)')' lm_size = ',paw_an(iatom)%lm_size 839 call wrtout(my_unt,msg,my_mode) 840 write(msg,'(a,i4)')' mesh_size = ',paw_an(iatom)%mesh_size 841 call wrtout(my_unt,msg,my_mode) 842 write(msg,'(a,i4)')' nkxc1 = ',paw_an(iatom)%nkxc1 843 call wrtout(my_unt,msg,my_mode) 844 write(msg,'(a,i4)')' nk3xc1 = ',paw_an(iatom)%nk3xc1 845 call wrtout(my_unt,msg,my_mode) 846 write(msg,'(a,i4)')' nspden = ',paw_an(iatom)%nspden 847 call wrtout(my_unt,msg,my_mode) 848 end do 849 850 end subroutine paw_an_print
m_paw_an/paw_an_redistribute [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_redistribute
FUNCTION
Redistribute an array of paw_an datastructures Input paw_an is given on a MPI communicator Output paw_an 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_an treated by current proc if not present, will be calculated in the present routine mpi_atmtab_out= --optional-- indexes of the output paw_an 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_an_out(:)]<type(paw_an_type)>= --optional-- if present, the redistributed datastructure does not replace the input one but is delivered in paw_an_out if not present, input and output datastructure are the same.
SIDE EFFECTS
paw_an(:)<type(paw_an_type)>= input (and eventually output) paw_an datastructures
PARENTS
m_paral_pert
CHILDREN
SOURCE
1445 subroutine paw_an_redistribute(paw_an,mpi_comm_in,mpi_comm_out,& 1446 & natom,mpi_atmtab_in,mpi_atmtab_out,paw_an_out,& 1447 & SendAtomProc,SendAtomList,RecvAtomProc,RecvAtomList) 1448 1449 1450 !This section has been created automatically by the script Abilint (TD). 1451 !Do not modify the following lines by hand. 1452 #undef ABI_FUNC 1453 #define ABI_FUNC 'paw_an_redistribute' 1454 !End of the abilint section 1455 1456 implicit none 1457 1458 !Arguments ------------------------------------ 1459 !scalars 1460 integer,intent(in) :: mpi_comm_in,mpi_comm_out 1461 integer,optional,intent(in) :: natom 1462 !arrays 1463 integer,intent(in),optional,target :: mpi_atmtab_in(:),mpi_atmtab_out(:) 1464 type(paw_an_type),allocatable,intent(inout) :: paw_an(:) 1465 type(paw_an_type),pointer,optional :: paw_an_out(:) !vz_i 1466 integer,intent(in),optional :: SendAtomProc(:),SendAtomList(:),RecvAtomProc(:),RecvAtomList(:) 1467 1468 !Local variables------------------------------- 1469 !scalars 1470 1471 integer :: algo_option,i1,iat_in,iat_out,iatom,ierr,iircv,iisend,imsg,imsg_current,imsg1 1472 integer :: iproc_rcv,iproc_send,ireq,me_exch,mpi_comm_exch,my_natom_in,my_natom_out,my_tag,natom_tot,nb_msg 1473 integer :: nb_dp,nb_int,nbmsg_incoming,nbrecvmsg,nbsend,nbsendreq,nbsent,nbrecv,next,npaw_an_sent 1474 integer :: nproc_in,nproc_out 1475 logical :: flag,in_place,message_yet_prepared,my_atmtab_in_allocated,my_atmtab_out_allocated,paral_atom 1476 !arrays 1477 integer :: buf_size(3),request1(3) 1478 integer,pointer :: my_atmtab_in(:),my_atmtab_out(:) 1479 integer,allocatable :: atmtab_send(:),atm_indx_in(:),atm_indx_out(:),buf_int1(:),From(:),request(:) 1480 integer,allocatable,target:: buf_int(:) 1481 integer,pointer :: buf_ints(:) 1482 logical, allocatable :: msg_pick(:) 1483 real(dp),allocatable :: buf_dp1(:) 1484 real(dp),allocatable,target :: buf_dp(:) 1485 real(dp),pointer :: buf_dps(:) 1486 type(coeffi1_type),target,allocatable :: tab_buf_int(:),tab_buf_atom(:) 1487 type(coeff1_type),target,allocatable :: tab_buf_dp(:) 1488 type(paw_an_type),allocatable :: paw_an_all(:) 1489 type(paw_an_type),pointer :: paw_an_out1(:) 1490 1491 ! ************************************************************************* 1492 1493 !@paw_an_type 1494 1495 in_place=(.not.present(paw_an_out)) 1496 my_natom_in=size(paw_an) 1497 1498 !If not "in_place", destroy the output datastructure 1499 if (.not.in_place) then 1500 if (associated(paw_an_out)) then 1501 call paw_an_free(paw_an_out) 1502 LIBPAW_DATATYPE_DEALLOCATE(paw_an_out) 1503 end if 1504 end if 1505 1506 !Special sequential case 1507 if (mpi_comm_in==xmpi_comm_self.and.mpi_comm_out==xmpi_comm_self) then 1508 if ((.not.in_place).and.(my_natom_in>0)) then 1509 LIBPAW_DATATYPE_ALLOCATE(paw_an_out,(my_natom_in)) 1510 call paw_an_nullify(paw_an_out) 1511 call paw_an_copy(paw_an,paw_an_out) 1512 end if 1513 return 1514 end if 1515 1516 !Get total natom 1517 if (present(natom)) then 1518 natom_tot=natom 1519 else 1520 natom_tot=my_natom_in 1521 call xmpi_sum(natom_tot,mpi_comm_in,ierr) 1522 end if 1523 1524 !Select input distribution 1525 if (present(mpi_atmtab_in)) then 1526 my_atmtab_in => mpi_atmtab_in 1527 my_atmtab_in_allocated=.false. 1528 else 1529 call get_my_atmtab(mpi_comm_in,my_atmtab_in,my_atmtab_in_allocated,& 1530 & paral_atom,natom_tot,my_natom_in) 1531 end if 1532 1533 !Select output distribution 1534 if (present(mpi_atmtab_out)) then 1535 my_natom_out=size(mpi_atmtab_out) 1536 my_atmtab_out => mpi_atmtab_out 1537 my_atmtab_out_allocated=.false. 1538 else 1539 call get_my_natom(mpi_comm_out,my_natom_out,natom_tot) 1540 call get_my_atmtab(mpi_comm_out,my_atmtab_out,my_atmtab_out_allocated,& 1541 & paral_atom,natom_tot) 1542 end if 1543 1544 !Select algo according to optional input arguments 1545 algo_option=1 1546 if (present(SendAtomProc).and.present(SendAtomList).and.& 1547 & present(RecvAtomProc).and.present(RecvAtomList)) algo_option=2 1548 1549 1550 !Brute force algorithm (allgather + scatter) 1551 !--------------------------------------------------------- 1552 if (algo_option==1) then 1553 1554 LIBPAW_DATATYPE_ALLOCATE(paw_an_all,(natom_tot)) 1555 call paw_an_nullify(paw_an_all) 1556 call paw_an_copy(paw_an,paw_an_all,comm_atom=mpi_comm_in,mpi_atmtab=my_atmtab_in) 1557 if (in_place) then 1558 call paw_an_free(paw_an) 1559 LIBPAW_DATATYPE_DEALLOCATE(paw_an) 1560 LIBPAW_DATATYPE_ALLOCATE(paw_an,(my_natom_out)) 1561 call paw_an_nullify(paw_an) 1562 call paw_an_copy(paw_an_all,paw_an,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 1563 else 1564 LIBPAW_DATATYPE_ALLOCATE(paw_an_out,(my_natom_out)) 1565 call paw_an_nullify(paw_an_out) 1566 call paw_an_copy(paw_an_all,paw_an_out,comm_atom=mpi_comm_out,mpi_atmtab=my_atmtab_out) 1567 end if 1568 call paw_an_free(paw_an_all) 1569 LIBPAW_DATATYPE_DEALLOCATE(paw_an_all) 1570 1571 1572 !Asynchronous algorithm (asynchronous communications) 1573 !--------------------------------------------------------- 1574 else if (algo_option==2) then 1575 1576 nbsend=size(SendAtomProc) ; nbrecv=size(RecvAtomProc) 1577 1578 if (in_place) then 1579 if (my_natom_out > 0) then 1580 LIBPAW_DATATYPE_ALLOCATE(paw_an_out1,(my_natom_out)) 1581 call paw_an_nullify(paw_an_out1) 1582 else 1583 LIBPAW_DATATYPE_ALLOCATE(paw_an_out1,(0)) 1584 end if 1585 else 1586 LIBPAW_DATATYPE_ALLOCATE(paw_an_out,(my_natom_out)) 1587 call paw_an_nullify(paw_an_out) 1588 paw_an_out1=>paw_an_out 1589 end if 1590 1591 nproc_in=xmpi_comm_size(mpi_comm_in) 1592 nproc_out=xmpi_comm_size(mpi_comm_out) 1593 if (nproc_in<=nproc_out) mpi_comm_exch=mpi_comm_out 1594 if (nproc_in>nproc_out) mpi_comm_exch=mpi_comm_in 1595 me_exch=xmpi_comm_rank(mpi_comm_exch) 1596 1597 ! Dimension put to the maximum to send 1598 LIBPAW_ALLOCATE(atmtab_send,(nbsend)) 1599 LIBPAW_ALLOCATE(atm_indx_in,(natom_tot)) 1600 atm_indx_in=-1 1601 do iatom=1,my_natom_in 1602 atm_indx_in(my_atmtab_in(iatom))=iatom 1603 end do 1604 LIBPAW_ALLOCATE(atm_indx_out,(natom_tot)) 1605 atm_indx_out=-1 1606 do iatom=1,my_natom_out 1607 atm_indx_out(my_atmtab_out(iatom))=iatom 1608 end do 1609 1610 LIBPAW_DATATYPE_ALLOCATE(tab_buf_int,(nbsend)) 1611 LIBPAW_DATATYPE_ALLOCATE(tab_buf_dp,(nbsend)) 1612 LIBPAW_DATATYPE_ALLOCATE(tab_buf_atom,(nbsend)) 1613 LIBPAW_ALLOCATE(request,(3*nbsend)) 1614 1615 ! A send buffer in an asynchrone communication couldn't be deallocate before it has been receive 1616 nbsent=0 ; ireq=0 ; iisend=0 ; nbsendreq=0 ; nb_msg=0 1617 do iisend=1,nbsend 1618 iproc_rcv=SendAtomProc(iisend) 1619 next=-1 1620 if (iisend < nbsend) next=SendAtomProc(iisend+1) 1621 if (iproc_rcv /= me_exch) then 1622 nbsent=nbsent+1 1623 atmtab_send(nbsent)=SendAtomList(iisend) ! we groups the atoms sends to the same process 1624 if (iproc_rcv /= next) then 1625 if (nbsent > 0) then 1626 ! Check if message has been yet prepared 1627 message_yet_prepared=.false. 1628 do imsg=1,nb_msg 1629 if (size(tab_buf_atom(imsg)%value) /= nbsent) then 1630 cycle 1631 else 1632 do imsg1=1,nbsent 1633 if (tab_buf_atom(imsg)%value(imsg1)/=atmtab_send(imsg1)) exit 1634 message_yet_prepared=.true. 1635 imsg_current=imsg 1636 end do 1637 end if 1638 end do 1639 ! Create the message 1640 if (.not.message_yet_prepared) then 1641 nb_msg=nb_msg+1 1642 call paw_an_isendreceive_fillbuffer( & 1643 & paw_an,atmtab_send,atm_indx_in,nbsent,buf_int,nb_int,buf_dp,nb_dp) 1644 LIBPAW_ALLOCATE(tab_buf_int(nb_msg)%value,(nb_int)) 1645 LIBPAW_ALLOCATE(tab_buf_dp(nb_msg)%value,(nb_dp)) 1646 tab_buf_int(nb_msg)%value(1:nb_int)=buf_int(1:nb_int) 1647 tab_buf_dp(nb_msg)%value(1:nb_dp)=buf_dp(1:nb_dp) 1648 LIBPAW_DEALLOCATE(buf_int) 1649 LIBPAW_DEALLOCATE(buf_dp) 1650 LIBPAW_ALLOCATE(tab_buf_atom(nb_msg)%value, (nbsent)) 1651 tab_buf_atom(nb_msg)%value(1:nbsent)=atmtab_send(1:nbsent) 1652 imsg_current=nb_msg 1653 end if 1654 ! Communicate 1655 buf_size(1)=size(tab_buf_int(imsg_current)%value) 1656 buf_size(2)=size(tab_buf_dp(imsg_current)%value) 1657 buf_size(3)=nbsent 1658 buf_ints=>tab_buf_int(imsg_current)%value 1659 buf_dps=>tab_buf_dp(imsg_current)%value 1660 my_tag=300 1661 ireq=ireq+1 1662 call xmpi_isend(buf_size,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 1663 my_tag=301 1664 ireq=ireq+1 1665 call xmpi_isend(buf_ints,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 1666 my_tag=302 1667 ireq=ireq+1 1668 call xmpi_isend(buf_dps,iproc_rcv,my_tag,mpi_comm_exch,request(ireq),ierr) 1669 nbsendreq=ireq 1670 nbsent=0 1671 end if 1672 end if 1673 else ! Just a renumbering, not a sending 1674 iat_in=atm_indx_in(SendAtomList(iisend)) 1675 iat_out=atm_indx_out(my_atmtab_in(iat_in)) 1676 call paw_an_copy(paw_an(iat_in:iat_in),paw_an_out1(iat_out:iat_out)) 1677 nbsent=0 1678 end if 1679 end do 1680 1681 LIBPAW_ALLOCATE(From,(nbrecv)) 1682 From(:)=-1 ; nbrecvmsg=0 1683 do iircv=1,nbrecv 1684 iproc_send=RecvAtomProc(iircv) !receive from ( RcvAtomProc is sorted by growing process ) 1685 next=-1 1686 if (iircv < nbrecv) next=RecvAtomProc(iircv+1) 1687 if (iproc_send /= me_exch .and. iproc_send/=next) then 1688 nbrecvmsg=nbrecvmsg+1 1689 From(nbrecvmsg)=iproc_send 1690 end if 1691 end do 1692 1693 LIBPAW_ALLOCATE(msg_pick,(nbrecvmsg)) 1694 msg_pick=.false. 1695 nbmsg_incoming=nbrecvmsg 1696 do while (nbmsg_incoming > 0) 1697 do i1=1,nbrecvmsg 1698 if (.not.msg_pick(i1)) then 1699 iproc_send=From(i1) 1700 flag=.false. 1701 my_tag=300 1702 call xmpi_iprobe(iproc_send,my_tag,mpi_comm_exch,flag,ierr) 1703 if (flag) then 1704 msg_pick(i1)=.true. 1705 call xmpi_irecv(buf_size,iproc_send,my_tag,mpi_comm_exch,request1(1),ierr) 1706 call xmpi_wait(request1(1),ierr) 1707 nb_int=buf_size(1) 1708 nb_dp=buf_size(2) 1709 npaw_an_sent=buf_size(3) 1710 LIBPAW_ALLOCATE(buf_int1,(nb_int)) 1711 LIBPAW_ALLOCATE(buf_dp1 ,(nb_dp)) 1712 my_tag=301 1713 call xmpi_irecv(buf_int1,iproc_send,my_tag,mpi_comm_exch,request1(2),ierr) 1714 my_tag=302 1715 call xmpi_irecv(buf_dp1,iproc_send,my_tag,mpi_comm_exch,request1(3),ierr) 1716 call xmpi_waitall(request1(2:3),ierr) 1717 call paw_an_isendreceive_getbuffer(paw_an_out1,npaw_an_sent,atm_indx_out,buf_int1,buf_dp1) 1718 nbmsg_incoming=nbmsg_incoming-1 1719 LIBPAW_DEALLOCATE(buf_int1) 1720 LIBPAW_DEALLOCATE(buf_dp1) 1721 end if 1722 end if 1723 end do 1724 end do 1725 LIBPAW_DEALLOCATE(msg_pick) 1726 1727 if (in_place) then 1728 call paw_an_free(paw_an) 1729 LIBPAW_DATATYPE_DEALLOCATE(paw_an) 1730 LIBPAW_DATATYPE_ALLOCATE(paw_an,(my_natom_out)) 1731 call paw_an_nullify(paw_an) 1732 call paw_an_copy(paw_an_out1,paw_an) 1733 call paw_an_free(paw_an_out1) 1734 LIBPAW_DATATYPE_DEALLOCATE(paw_an_out1) 1735 end if 1736 1737 ! Wait for deallocating arrays that all sending operations has been realized 1738 if (nbsendreq > 0) then 1739 call xmpi_waitall(request(1:nbsendreq),ierr) 1740 end if 1741 1742 ! Deallocate buffers 1743 do i1=1,nb_msg 1744 LIBPAW_DEALLOCATE(tab_buf_int(i1)%value) 1745 LIBPAW_DEALLOCATE(tab_buf_dp(i1)%value) 1746 LIBPAW_DEALLOCATE(tab_buf_atom(i1)%value) 1747 end do 1748 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_int) 1749 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_dp) 1750 LIBPAW_DATATYPE_DEALLOCATE(tab_buf_atom) 1751 LIBPAW_DEALLOCATE(From) 1752 LIBPAW_DEALLOCATE(request) 1753 LIBPAW_DEALLOCATE(atmtab_send) 1754 LIBPAW_DEALLOCATE(atm_indx_in) 1755 LIBPAW_DEALLOCATE(atm_indx_out) 1756 1757 end if !algo_option 1758 1759 !Eventually release temporary pointers 1760 call free_my_atmtab(my_atmtab_in,my_atmtab_in_allocated) 1761 call free_my_atmtab(my_atmtab_out,my_atmtab_out_allocated) 1762 1763 end subroutine paw_an_redistribute
m_paw_an/paw_an_reset_flags [ Functions ]
[ Top ] [ m_paw_an ] [ Functions ]
NAME
paw_an_reset_flags
FUNCTION
Set all paw_an flags to 1 (force the recomputation of all arrays)
SIDE EFFECTS
Paw_an<type(Paw_an_type)>=paw_an datastructure
PARENTS
dfpt_nstpaw,dfpt_scfcv,scfcv
CHILDREN
SOURCE
1785 subroutine paw_an_reset_flags(Paw_an) 1786 1787 1788 !This section has been created automatically by the script Abilint (TD). 1789 !Do not modify the following lines by hand. 1790 #undef ABI_FUNC 1791 #define ABI_FUNC 'paw_an_reset_flags' 1792 !End of the abilint section 1793 1794 implicit none 1795 1796 !Arguments ------------------------------------ 1797 !arrays 1798 type(Paw_an_type),intent(inout) :: Paw_an(:) 1799 1800 !Local variables------------------------------- 1801 integer :: iat,natom 1802 1803 ! ************************************************************************* 1804 1805 !@Paw_an_type 1806 1807 natom=SIZE(Paw_an);if (natom==0) return 1808 do iat=1,natom 1809 if (Paw_an(iat)%has_kxc >0) Paw_an(iat)%has_kxc =1 1810 if (Paw_an(iat)%has_k3xc >0) Paw_an(iat)%has_k3xc =1 1811 if (Paw_an(iat)%has_vhartree>0) Paw_an(iat)%has_vhartree=1 1812 if (Paw_an(iat)%has_vxc >0) Paw_an(iat)%has_vxc =1 1813 if (Paw_an(iat)%has_vxcval >0) Paw_an(iat)%has_vxcval =1 1814 if (Paw_an(iat)%has_vxc_ex >0) Paw_an(iat)%has_vxc_ex =1 1815 end do 1816 1817 end subroutine paw_an_reset_flags
m_paw_an/paw_an_type [ Types ]
[ Top ] [ m_paw_an ] [ Types ]
NAME
paw_an_type
FUNCTION
For PAW, various arrays given on ANgular mesh or ANgular moments
SOURCE
67 type,public :: paw_an_type 68 69 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 70 ! declared in another part of ABINIT, that might need to take into account your modification. 71 72 !Integer scalars 73 74 integer :: angl_size 75 ! Dimension of paw angular mesh (angl_size=ntheta*nphi) 76 77 integer :: cplex 78 ! cplex=1 if potentials/densities are real, 2 if they are complex 79 80 integer :: has_kxc 81 ! set to 1 if xc kernels kxc1 and kxct1 are allocated and used 82 ! 2 if they are already computed 83 84 integer :: has_k3xc 85 ! set to 1 if xc kernel derivatives k3xc1 and k3xct1 are allocated and used 86 ! 2 if it is already computed 87 88 integer :: has_vhartree 89 ! set to 1 if vh1 and vht1 are allocated and used 90 ! 2 if they are already computed 91 92 integer :: has_vxc 93 ! set to 1 if vxc1 and vxct1 are allocated and used 94 ! 2 if they are already computed 95 96 integer :: has_vxcval 97 ! set to 1 if vxc1_val and vxct1_val are allocated and used 98 ! 2 if they are already computed 99 100 integer :: has_vxc_ex 101 ! set to 1 if vxc_ex and is allocated and used 102 ! 2 if it is already computed 103 104 integer :: itypat 105 ! itypat=type of the atom 106 107 integer :: lm_size 108 ! lm_size=(l_size)**2 109 ! l is Maximum value of l+1 leading to non zero Gaunt coeffs (l_size=2*l_max+1) 110 111 integer :: mesh_size 112 ! Dimension of radial mesh for arrays contained in this paw_an datastructure 113 ! May be different from pawrad%mesh_size 114 115 integer :: nkxc1 116 ! number of independent components of Kxc1 and Kxct1 117 ! (usually 3 for LDA, 23 for GGA) 118 119 integer :: nk3xc1 120 ! number of independent components of K3xc1 and K3xct1 121 ! (usually 4 for LDA, not available for GGA) 122 123 integer :: nspden 124 ! Number of spin-density components 125 126 !Logical arrays 127 128 logical, allocatable :: lmselect(:) 129 ! lmselect(lm_size) 130 ! lmselect(ilm)=select the non-zero LM-moments of "one-center" densities/potentials 131 132 !Real (real(dp)) arrays 133 134 real(dp), allocatable :: kxc1 (:,:,:) 135 ! kxc1(cplex*mesh_size,lm_size or angl_size,nkxc1) 136 ! Gives xc kernel inside the sphere 137 ! (theta,phi) values of kernel if pawxcdev=0 138 ! LM-moments of kernel if pawxcdev/=0 139 140 real(dp), allocatable :: kxct1 (:,:,:) 141 ! kxct1(cplex*mesh_size,lm_size or angl_size,nkxc1) 142 ! Gives xc pseudo kernel inside the sphere 143 ! (theta,phi) values of kernel if pawxcdev=0 144 ! LM-moments of kernel if pawxcdev/=0 145 146 real(dp), allocatable :: k3xc1 (:,:,:) 147 ! k3xc1(cplex*mesh_size,lm_size or angl_size,nk3xc1) 148 ! Gives xc kernel derivative inside the sphere 149 ! (theta,phi) values of kernel derivative if pawxcdev=0 150 ! LM-moments of kernel derivative if pawxcdev/=0 => NOT AVAILABLE YET 151 152 real(dp), allocatable :: k3xct1 (:,:,:) 153 ! k3xct1(cplex*mesh_size,lm_size or angl_size,nk3xc1) 154 ! Gives xc pseudo kernel derivative inside the sphere 155 ! (theta,phi) values of kernel derivative if pawxcdev=0 156 ! LM-moments of kernel derivative if pawxcdev/=0 => NOT AVAILABLE YET 157 158 real(dp), allocatable :: vh1 (:,:,:) 159 ! vh1(cplex*mesh_size,lm_size,nspden) 160 ! Gives Hartree potential LM-moments inside the sphere 161 162 real(dp), allocatable :: vht1 (:,:,:) 163 ! vht1(cplex*mesh_size,lm_size,nspden) 164 ! Gives Hartree tilde potential LM-moments inside the sphere 165 166 real(dp), allocatable :: vxc1 (:,:,:) 167 ! vxc1(cplex*mesh_size,lm_size or angl_size,nspden) 168 ! Gives xc potential inside the sphere 169 ! (theta,phi) values of potential if pawxcdev=0 170 ! LM-moments of potential if pawxcdev/=0 171 172 real(dp), allocatable :: vxc1_val (:,:,:) 173 ! vxc1_val(cplex*mesh_size,lm_size or angl_size,nspden) (Usually real, Mainly used for GW) 174 ! Gives xc potential inside the sphere arising from valence only electrons 175 ! (theta,phi) values of potential if pawxcdev=0 176 ! LM-moments of potential if pawxcdev/=0 177 178 real(dp), allocatable :: vxct1 (:,:,:) 179 ! vxct1(cplex*mesh_size,angl_size,nspden) 180 ! Gives xc pseudo potential inside the sphere 181 ! (theta,phi) values of potential if pawxcdev=0 182 ! LM-moments of potential if pawxcdev/=0 183 184 real(dp), allocatable :: vxct1_val (:,:,:) 185 ! vxct1_val(cplex*mesh_size,angl_size,nspden) (Usually real, Mainly used for GW) 186 ! Gives xc pseudo potential inside the sphere 187 ! (theta,phi) values of potential if pawxcdev=0 188 ! LM-moments of potential if pawxcdev/=0 189 190 real(dp), allocatable :: vxc_ex (:,:,:) 191 ! vxc_ex(cplex*mesh_size,angl_size,nspden) 192 ! Gives xc potential for local exact exchange inside the sphere 193 ! (theta,phi) values of potential if pawxcdev=0 194 ! LM-moments of potential if pawxcdev/=0 195 196 end type paw_an_type