TABLE OF CONTENTS
- abinit/m_pawcprj
- m_pawcprj/paw_overlap
- m_pawcprj/pawcprj_alloc
- m_pawcprj/pawcprj_axpby
- m_pawcprj/pawcprj_bcast
- m_pawcprj/pawcprj_conjg
- m_pawcprj/pawcprj_copy
- m_pawcprj/pawcprj_free
- m_pawcprj/pawcprj_gather_spin
- m_pawcprj/pawcprj_get
- m_pawcprj/pawcprj_getdim
- m_pawcprj/pawcprj_lincom
- m_pawcprj/pawcprj_mpi_allgather
- m_pawcprj/pawcprj_mpi_exch
- m_pawcprj/pawcprj_mpi_recv
- m_pawcprj/pawcprj_mpi_send
- m_pawcprj/pawcprj_mpi_sum
- m_pawcprj/pawcprj_output
- m_pawcprj/pawcprj_put
- m_pawcprj/pawcprj_reorder
- m_pawcprj/pawcprj_set_zero
- m_pawcprj/pawcprj_symkn
- m_pawcprj/pawcprj_transpose
- m_pawcprj/pawcprj_type
- m_pawcprj/pawcprj_zaxpby
abinit/m_pawcprj [ Modules ]
[ Top ] [ abinit ] [ Modules ]
NAME
m_pawcprj
FUNCTION
This module contains functions used to manipulate variables of structured datatype pawcprj_type. pawcprj_type variables are <p_lmn|Cnk> projected quantities, where |p_lmn> are non-local projectors |Cnk> are wave functions
COPYRIGHT
Copyright (C) 2012-2018 ABINIT group (MT,JWZ) 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
25 #include "libpaw.h" 26 27 module m_pawcprj 28 29 USE_DEFS 30 USE_MSG_HANDLING 31 USE_MPI_WRAPPERS 32 USE_MEMORY_PROFILING 33 34 use m_pawtab, only : pawtab_type 35 36 implicit none 37 38 private
m_pawcprj/paw_overlap [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
paw_overlap
FUNCTION
Helper function returning the onsite contribution to the overlap between two states.
INPUTS
spinor_comm= (optional) communicator over spinorial components typat(:)=The type of each atom. Pawtab(ntypat)<type(pawtab_type)>=paw tabulated starting data. cprj1,cprj2<pawcprj_type> Projected wave functions <Proj_i|Cnk> with all NL projectors for the left and the right wavefunction,respectively.
OUTPUT
PARENTS
CHILDREN
xmpi_sum
SOURCE
2965 function paw_overlap(cprj1,cprj2,typat,pawtab,spinor_comm) result(onsite) 2966 2967 2968 !This section has been created automatically by the script Abilint (TD). 2969 !Do not modify the following lines by hand. 2970 #undef ABI_FUNC 2971 #define ABI_FUNC 'paw_overlap' 2972 !End of the abilint section 2973 2974 implicit none 2975 2976 !Arguments ------------------------------------ 2977 !scalars 2978 integer,intent(in),optional :: spinor_comm 2979 !arrays 2980 integer,intent(in) :: typat(:) 2981 real(dp) :: onsite(2) 2982 type(pawcprj_type),intent(in) :: cprj1(:,:),cprj2(:,:) 2983 type(pawtab_type),intent(in) :: pawtab(:) 2984 2985 !Local variables------------------------------- 2986 !scalars 2987 integer :: iatom,ilmn,itypat,j0lmn,jlmn,klmn,natom,nspinor,isp 2988 real(dp) :: sij 2989 character(len=500) :: msg 2990 !arrays 2991 2992 ! ************************************************************************* 2993 2994 natom=SIZE(typat) 2995 2996 if (SIZE(cprj1,DIM=1)/=SIZE(cprj2,DIM=1) .or. SIZE(cprj1,DIM=1)/=natom) then 2997 write(msg,'(a,3i4)')' Wrong size in typat, cprj1, cprj2 : ',natom,SIZE(cprj1),SIZE(cprj2) 2998 MSG_ERROR(msg) 2999 end if 3000 3001 nspinor = SIZE(cprj1,DIM=2) 3002 3003 onsite=zero 3004 do iatom=1,natom 3005 itypat=typat(iatom) 3006 do jlmn=1,pawtab(itypat)%lmn_size 3007 j0lmn=jlmn*(jlmn-1)/2 3008 do ilmn=1,jlmn 3009 klmn=j0lmn+ilmn 3010 sij=pawtab(itypat)%sij(klmn); if (jlmn==ilmn) sij=sij*half 3011 if (ABS(sij)>tol16) then 3012 do isp=1,nspinor 3013 3014 onsite(1)=onsite(1) + sij*( & 3015 & cprj1(iatom,isp)%cp(1,ilmn) * cprj2(iatom,isp)%cp(1,jlmn) & 3016 & +cprj1(iatom,isp)%cp(2,ilmn) * cprj2(iatom,isp)%cp(2,jlmn) & 3017 & +cprj1(iatom,isp)%cp(1,jlmn) * cprj2(iatom,isp)%cp(1,ilmn) & 3018 & +cprj1(iatom,isp)%cp(2,jlmn) * cprj2(iatom,isp)%cp(2,ilmn) & 3019 & ) 3020 3021 onsite(2)=onsite(2) + sij*( & 3022 & cprj1(iatom,isp)%cp(1,ilmn) * cprj2(iatom,isp)%cp(2,jlmn) & 3023 & -cprj1(iatom,isp)%cp(2,ilmn) * cprj2(iatom,isp)%cp(1,jlmn) & 3024 & +cprj1(iatom,isp)%cp(1,jlmn) * cprj2(iatom,isp)%cp(2,ilmn) & 3025 & -cprj1(iatom,isp)%cp(2,jlmn) * cprj2(iatom,isp)%cp(1,ilmn) & 3026 & ) 3027 end do 3028 end if 3029 end do 3030 end do 3031 end do 3032 3033 if (present(spinor_comm)) then 3034 call xmpi_sum(onsite,spinor_comm,isp) 3035 end if 3036 3037 end function paw_overlap
m_pawcprj/pawcprj_alloc [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_alloc
FUNCTION
Allocation of a cprj datastructure
INPUTS
ncpgr=number of gradients to be allocated nlmn(:)=sizes of cprj%cp
SIDE EFFECTS
cprj(:,:) <type(pawcprj_type)>= cprj datastructure
PARENTS
berryphase_new,calc_optical_mels,calc_sigc_me,calc_sigx_me,calc_vhxc_me calc_wf_qp,cchi0,cchi0q0,cchi0q0_intraband,cgwf,chebfi,chern_number classify_bands,cohsex_me,ctocprj,d2frnl,m_datafordmft,debug_tools dfpt_accrho,dfpt_cgwf,dfpt_looppert,dfpt_nstpaw,dfpt_scfcv,dfpt_vtowfk dfpt_wfkfermi,dotprod_set_cgcprj,dotprodm_sumdiag_cgcprj,energy exc_build_block,exc_build_ham,exc_plot,extrapwf,fock2ACE,forstr forstrnps,getgh1c,getgh2c,getghc,getgsc,initberry,ks_ddiago lincom_cgcprj,m_electronpositron,m_fock,m_invovl,m_io_kss,m_pawcprj m_plowannier,m_shirley,m_wfd,make_grad_berry,nonlop,optics_paw optics_paw_core,outkss,partial_dos_fractions_paw,paw_symcprj,pawmkaewf pawmkrhoij,posdoppler,prep_calc_ucrpa,rf2_init,scfcv,setup_positron sigma,smatrix_pawinit,suscep_stat,update_e_field_vars,vtorho,vtowfk wf_mixing,wfd_pawrhoij,wfd_vnlpsi,wvl_hpsitopsi
CHILDREN
xmpi_sum
SOURCE
140 subroutine pawcprj_alloc(cprj,ncpgr,nlmn) 141 142 143 !This section has been created automatically by the script Abilint (TD). 144 !Do not modify the following lines by hand. 145 #undef ABI_FUNC 146 #define ABI_FUNC 'pawcprj_alloc' 147 !End of the abilint section 148 149 implicit none 150 151 !Arguments ------------------------------------ 152 !scalars 153 integer,intent(in) :: ncpgr 154 !arrays 155 integer,intent(in) :: nlmn(:) 156 type(pawcprj_type),intent(inout) :: cprj(:,:) 157 158 !Local variables------------------------------- 159 !scalars 160 integer :: ii,jj,n1dim,n2dim,nn 161 character(len=500) :: msg 162 163 ! ************************************************************************* 164 165 n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2);nn=size(nlmn,dim=1) 166 if (nn/=n1dim) then 167 write(msg,*) 'wrong sizes (pawcprj_alloc)! :',nn,n1dim 168 MSG_ERROR(msg) 169 end if 170 171 do jj=1,n2dim 172 do ii=1,n1dim 173 if (allocated(cprj(ii,jj)%cp)) then 174 LIBPAW_DEALLOCATE(cprj(ii,jj)%cp) 175 end if 176 if (allocated(cprj(ii,jj)%dcp)) then 177 LIBPAW_DEALLOCATE(cprj(ii,jj)%dcp) 178 end if 179 nn=nlmn(ii) 180 cprj(ii,jj)%nlmn=nn 181 LIBPAW_ALLOCATE(cprj(ii,jj)%cp,(2,nn)) 182 cprj(ii,jj)%cp=zero 183 cprj(ii,jj)%ncpgr=ncpgr 184 if (ncpgr>0) then 185 LIBPAW_ALLOCATE(cprj(ii,jj)%dcp,(2,ncpgr,nn)) 186 cprj(ii,jj)%dcp=zero 187 end if 188 end do 189 end do 190 191 end subroutine pawcprj_alloc
m_pawcprj/pawcprj_axpby [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_axpby
FUNCTION
Apply AXPBY (blas-like) operation with 2 cprj datastructures: cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:) alpha and beta are REAL scalars
INPUTS
alpha,beta= alpha,beta REAL factors cprjx(:,:) <type(pawcprj_type)>= input cprjx datastructure
SIDE EFFECTS
cprjy(:,:) <type(pawcprj_type)>= input/output cprjy datastructure
PARENTS
chebfi,dfpt_cgwf,dfpt_wfkfermi,getdc1,m_invovl,wf_mixing
CHILDREN
xmpi_sum
SOURCE
470 subroutine pawcprj_axpby(alpha,beta,cprjx,cprjy) 471 472 473 !This section has been created automatically by the script Abilint (TD). 474 !Do not modify the following lines by hand. 475 #undef ABI_FUNC 476 #define ABI_FUNC 'pawcprj_axpby' 477 !End of the abilint section 478 479 implicit none 480 481 !Arguments ------------------------------------ 482 !scalars 483 real(dp),intent(in) :: alpha,beta 484 !arrays 485 type(pawcprj_type),intent(in) :: cprjx(:,:) 486 type(pawcprj_type),intent(inout) :: cprjy(:,:) 487 488 !Local variables------------------------------- 489 !scalars 490 integer :: ii,jj,kk,n1dimx,n1dimy,n2dimx,n2dimy,ncpgrx,ncpgry,nlmn 491 character(len=500) :: msg 492 493 ! ************************************************************************* 494 495 n1dimy=size(cprjy,dim=1);n2dimy=size(cprjy,dim=2);ncpgry=cprjy(1,1)%ncpgr 496 if (abs(alpha)>tol16) then 497 n1dimx=size(cprjx,dim=1);n2dimx=size(cprjx,dim=2);ncpgrx=cprjx(1,1)%ncpgr 498 msg = "" 499 if (n1dimx/=n1dimy) msg = TRIM(msg)//"Error in pawcprj_axpby: n1 wrong sizes !"//ch10 500 if (n2dimx/=n2dimy) msg = TRIM(msg)//"Error in pawcprj_axpby: n2 wrong sizes !"//ch10 501 if (ncpgrx/=ncpgry) msg = TRIM(msg)//"Error in pawcprj_axpby: ncpgr wrong sizes !"//ch10 502 if (LEN_TRIM(msg) > 0) then 503 MSG_ERROR(msg) 504 end if 505 else 506 n1dimx=0;n2dimx=0;ncpgrx=0 507 end if 508 509 if (abs(alpha)<=tol16) then 510 do jj=1,n2dimy 511 do ii=1,n1dimy 512 nlmn=cprjy(ii,jj)%nlmn 513 do kk=1,nlmn 514 cprjy(ii,jj)%cp(1:2,kk)=beta*cprjy(ii,jj)%cp(1:2,kk) 515 end do 516 end do 517 end do 518 if (ncpgry>0) then 519 do jj=1,n2dimy 520 do ii=1,n1dimy 521 nlmn=cprjy(ii,jj)%nlmn 522 do kk=1,nlmn 523 cprjy(ii,jj)%dcp(1:2,1:ncpgry,kk)=beta*cprjy(ii,jj)%dcp(1:2,1:ncpgry,kk) 524 end do 525 end do 526 end do 527 end if 528 else if (abs(beta)<=tol16) then 529 do jj=1,n2dimx 530 do ii=1,n1dimx 531 nlmn=cprjx(ii,jj)%nlmn 532 cprjy(ii,jj)%nlmn=nlmn 533 do kk=1,nlmn 534 cprjy(ii,jj)%cp(1:2,kk)=alpha*cprjx(ii,jj)%cp(1:2,kk) 535 end do 536 end do 537 end do 538 if (ncpgrx>0) then 539 do jj=1,n2dimx 540 do ii=1,n1dimx 541 nlmn=cprjx(ii,jj)%nlmn 542 do kk=1,nlmn 543 cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk)=alpha*cprjx(ii,jj)%dcp(1:2,1:ncpgrx,kk) 544 end do 545 end do 546 end do 547 end if 548 else ! alpha/=0 and beta/=0 549 do jj=1,n2dimx 550 do ii=1,n1dimx 551 nlmn=cprjx(ii,jj)%nlmn 552 cprjy(ii,jj)%nlmn=nlmn 553 do kk=1,nlmn 554 cprjy(ii,jj)%cp(1:2,kk)=alpha*cprjx(ii,jj)%cp(1:2,kk) & 555 & +beta *cprjy(ii,jj)%cp(1:2,kk) 556 end do 557 end do 558 end do 559 if (ncpgrx>0) then 560 do jj=1,n2dimx 561 do ii=1,n1dimx 562 nlmn=cprjx(ii,jj)%nlmn 563 do kk=1,nlmn 564 cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk)=alpha*cprjx(ii,jj)%dcp(1:2,1:ncpgrx,kk) & 565 & +beta *cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk) 566 end do 567 end do 568 end do 569 end if 570 end if 571 572 end subroutine pawcprj_axpby
m_pawcprj/pawcprj_bcast [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_bcast
FUNCTION
Broadcast a pawcprj_type from master to all nodes inside a MPI communicator.
INPUTS
natom=Number of atoms (size of the first dimension of Cprj). n2dim=Size of the second dimension of Cprj. ncpgr=Number of gradients that have to be cast. It is a bit redundant but, it can be used to broad cast only the %cp"s without caring about the gradients. Just set it to 0 but be careful! nlmn(natom)=Number of nlm partial waves for each atom. master=ID of the sending node in spaceComm. spaceComm=MPI Communicator.
OUTPUT
ierr=Error status. Cprj(natom,n2dim)<pawcprj_type>=The datatype to be transmitted by master and received by the others nodes.
PARENTS
m_fock,posdoppler
CHILDREN
xmpi_sum
SOURCE
2375 subroutine pawcprj_bcast(Cprj,natom,n2dim,nlmn,ncpgr,master,spaceComm,ierr) 2376 2377 2378 !This section has been created automatically by the script Abilint (TD). 2379 !Do not modify the following lines by hand. 2380 #undef ABI_FUNC 2381 #define ABI_FUNC 'pawcprj_bcast' 2382 !End of the abilint section 2383 2384 implicit none 2385 2386 !Arguments ------------------------------------ 2387 !scalars 2388 integer,intent(in) :: natom,n2dim,ncpgr,master,spaceComm 2389 integer,intent(out) :: ierr 2390 !arrays 2391 integer,intent(in) :: nlmn(natom) 2392 type(pawcprj_type),intent(inout) :: Cprj(natom,n2dim) 2393 2394 !Local variables------------------------------- 2395 !scalars 2396 integer :: iat,jj,n1dim,nn 2397 integer :: ntotcp,ipck,rank,nprocs 2398 character(len=100) :: msg 2399 !arrays 2400 real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:) 2401 2402 ! ************************************************************************* 2403 2404 ierr=0 2405 nprocs = xmpi_comm_size(spaceComm) 2406 if (nprocs==1) return 2407 2408 rank = xmpi_comm_rank(spaceComm) 2409 2410 nn=size(nlmn,dim=1) 2411 n1dim=size(Cprj,dim=1) 2412 if (nn/=n1dim) then 2413 msg='size mismatch in natom (pawcprj_bcast)!' 2414 MSG_BUG(msg) 2415 end if 2416 2417 ntotcp=n2dim*SUM(nlmn(:)) 2418 2419 LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp)) 2420 if (ncpgr/=0) then 2421 LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp)) 2422 end if 2423 2424 !=== Master packs Cprj === 2425 !Write a routine to pack/unpack? 2426 if (rank==master) then 2427 ipck=0 2428 do jj=1,n2dim 2429 do iat=1,natom 2430 nn=nlmn(iat) 2431 buffer_cp(:,ipck+1:ipck+nn)=Cprj(iat,jj)%cp(:,1:nn) 2432 if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=Cprj(iat,jj)%dcp(:,:,1:nn) 2433 ipck=ipck+nn 2434 end do 2435 end do 2436 end if 2437 2438 !=== Transmit data === 2439 call xmpi_bcast(buffer_cp,master,spaceComm,ierr) 2440 if (ncpgr/=0) then 2441 call xmpi_bcast(buffer_cpgr,master,spaceComm,ierr) 2442 end if 2443 2444 !=== UnPack the received buffer === 2445 if (rank/=master) then 2446 ipck=0 2447 do jj=1,n2dim 2448 do iat=1,natom 2449 nn=nlmn(iat) 2450 Cprj(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn) 2451 if (ncpgr/=0) Cprj(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn) 2452 ipck=ipck+nn 2453 end do 2454 end do 2455 end if 2456 2457 LIBPAW_DEALLOCATE(buffer_cp) 2458 if (ncpgr/=0) then 2459 LIBPAW_DEALLOCATE(buffer_cpgr) 2460 end if 2461 2462 end subroutine pawcprj_bcast
m_pawcprj/pawcprj_conjg [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_conjg
FUNCTION
conjugate a cprj datastructures: cprj(:,:) <- conjugate(cprj(:,:))
INPUTS
SIDE EFFECTS
cprj(:,:) <type(pawcprj_type)>= input/output cprj datastructure
PARENTS
CHILDREN
xmpi_sum
SOURCE
918 subroutine pawcprj_conjg(cprj) 919 920 921 !This section has been created automatically by the script Abilint (TD). 922 !Do not modify the following lines by hand. 923 #undef ABI_FUNC 924 #define ABI_FUNC 'pawcprj_conjg' 925 !End of the abilint section 926 927 implicit none 928 929 !Arguments ------------------------------------ 930 !scalars 931 !arrays 932 type(pawcprj_type),intent(inout) :: cprj(:,:) 933 934 !Local variables------------------------------- 935 !scalars 936 integer :: ii,jj,kk,n1dim,n2dim,ncpgr,nlmn 937 !arrays 938 939 ! ************************************************************************* 940 941 942 n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2);ncpgr=cprj(1,1)%ncpgr 943 944 do jj=1,n2dim 945 do ii=1,n1dim 946 nlmn=cprj(ii,jj)%nlmn 947 do kk=1,nlmn 948 cprj(ii,jj)%cp(2,kk)=-cprj(ii,jj)%cp(2,kk) 949 end do 950 end do 951 end do 952 if (ncpgr>0) then 953 do jj=1,n2dim 954 do ii=1,n1dim 955 nlmn=cprj(ii,jj)%nlmn 956 do kk=1,nlmn 957 cprj(ii,jj)%dcp(2,1:ncpgr,kk)=-cprj(ii,jj)%dcp(2,1:ncpgr,kk) 958 end do 959 end do 960 end do 961 end if 962 963 end subroutine pawcprj_conjg
m_pawcprj/pawcprj_copy [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_copy
FUNCTION
Copy a cprj datastructure into another
INPUTS
icpgr= (optional argument) if present, only component icpgr of input cprj gradient is copied into output cprj Not used if cprj(:,:)%ncpgr<icpgr -1 only copy cp cprj_in(:,:) <type(pawcprj_type)>= input cprj datastructure
OUTPUT
cprj_out(:,:) <type(pawcprj_type)>= output cprj datastructure
NOTES
MG: What about an option to report a pointer to cprj_in?
PARENTS
berryphase_new,calc_sigc_me,calc_sigx_me,cchi0q0,cchi0q0_intraband,cgwf chebfi,classify_bands,cohsex_me,corrmetalwf1,dfpt_looppert,dfpt_nstpaw dfpt_vtowfk,dfpt_wfkfermi,extrapwf,getgsc,m_electronpositron,m_fock m_pawcprj,m_wfd,make_grad_berry,nonlop,outkss,paw_symcprj,posdoppler prep_calc_ucrpa,setup_positron,vtowfk
CHILDREN
xmpi_sum
SOURCE
354 subroutine pawcprj_copy(cprj_in,cprj_out,& 355 & icpgr) ! optional argument 356 357 358 !This section has been created automatically by the script Abilint (TD). 359 !Do not modify the following lines by hand. 360 #undef ABI_FUNC 361 #define ABI_FUNC 'pawcprj_copy' 362 !End of the abilint section 363 364 implicit none 365 366 !Arguments ------------------------------------ 367 !scalars 368 integer,intent(in),optional :: icpgr 369 !arrays 370 type(pawcprj_type),intent(in) :: cprj_in(:,:) 371 type(pawcprj_type),intent(inout) :: cprj_out(:,:) 372 373 !Local variables------------------------------- 374 !scalars 375 integer :: ii,jj,kk,n1dim_in,n1dim_out,n2dim_in,n2dim_out,ncpgr_in,ncpgr_out,nlmn 376 logical :: has_icpgr,copy_dcp 377 character(len=500) :: msg 378 379 ! ************************************************************************* 380 381 n1dim_in=size(cprj_in,dim=1); n1dim_out=size(cprj_out,dim=1) 382 n2dim_in=size(cprj_in,dim=2); n2dim_out=size(cprj_out,dim=2) 383 ncpgr_in=cprj_in(1,1)%ncpgr; ncpgr_out=cprj_out(1,1)%ncpgr 384 385 if (n1dim_in/=n1dim_out) then 386 write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: n1 wrong sizes ",n1dim_in,n1dim_out 387 MSG_ERROR(msg) 388 end if 389 if (n2dim_in/=n2dim_out) then 390 write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: n2 wrong sizes ",n2dim_in,n2dim_out 391 MSG_ERROR(msg) 392 end if 393 if (ncpgr_in<ncpgr_out) then 394 write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: ncpgr wrong sizes ",ncpgr_in,ncpgr_out 395 MSG_ERROR(msg) 396 end if 397 398 !Check if icgr is present and if dcp have to be copy 399 has_icpgr=present(icpgr) 400 copy_dcp = .TRUE. 401 if(has_icpgr)then 402 copy_dcp = icpgr>=0 403 end if 404 405 do jj=1,n2dim_in 406 do ii=1,n1dim_in 407 nlmn=cprj_in(ii,jj)%nlmn 408 cprj_out(ii,jj)%nlmn =nlmn 409 do kk=1,nlmn 410 cprj_out(ii,jj)%cp(1:2,kk)=cprj_in(ii,jj)%cp(1:2,kk) 411 end do 412 end do 413 end do 414 415 if (ncpgr_in>0.and.copy_dcp) then 416 if (has_icpgr) has_icpgr=(ncpgr_out>0.and.icpgr>0.or.icpgr<=ncpgr_in) 417 418 if (has_icpgr) then 419 do jj=1,n2dim_in 420 do ii=1,n1dim_in 421 nlmn=cprj_in(ii,jj)%nlmn 422 do kk=1,nlmn 423 cprj_out(ii,jj)%dcp(1:2,1,kk)=cprj_in(ii,jj)%dcp(1:2,icpgr,kk) 424 end do 425 end do 426 end do 427 else 428 if (ncpgr_out>=ncpgr_in) then 429 do jj=1,n2dim_in 430 do ii=1,n1dim_in 431 nlmn=cprj_in(ii,jj)%nlmn 432 do kk=1,nlmn 433 cprj_out(ii,jj)%dcp(1:2,1:ncpgr_in,kk)=cprj_in(ii,jj)%dcp(1:2,1:ncpgr_in,kk) 434 end do 435 end do 436 end do 437 end if 438 end if 439 end if 440 441 end subroutine pawcprj_copy
m_pawcprj/pawcprj_free [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_free
FUNCTION
Deallocation of a cprj datastructure
SIDE EFFECTS
cprj(:,:) <type(pawcprj_type)>= cprj datastructure
PARENTS
berryphase_new,calc_optical_mels,calc_sigc_me,calc_sigx_me,calc_vhxc_me calc_wf_qp,cchi0,cchi0q0,cchi0q0_intraband,cgwf,chebfi,chern_number classify_bands,cohsex_me,ctocprj,d2frnl,datafordmft,debug_tools dfpt_accrho,dfpt_cgwf,dfpt_looppert,dfpt_nstpaw,dfpt_scfcv,dfpt_vtowfk dfpt_wfkfermi,dotprod_set_cgcprj,dotprodm_sumdiag_cgcprj,energy exc_build_block,exc_build_ham,exc_plot,extrapwf,fock2ACE,forstr forstrnps,getgh1c,getgh2c,getghc,getgsc,ks_ddiago,lincom_cgcprj m_efield,m_electronpositron,m_fock,m_gkk,m_invovl,m_io_kss,m_pawcprj m_phgamma,m_phpi,m_plowannier,m_scf_history,m_shirley,m_sigmaph,m_wfd make_grad_berry,nonlop,optics_paw,optics_paw_core,outkss partial_dos_fractions_paw,paw_symcprj,pawmkaewf,pawmkrhoij,posdoppler prep_calc_ucrpa,rf2_init,scfcv,setup_positron,sigma,smatrix_pawinit suscep_stat,update_e_field_vars,vtorho,vtowfk,wf_mixing,wfd_pawrhoij wfd_vnlpsi
CHILDREN
xmpi_sum
SOURCE
227 subroutine pawcprj_free(cprj) 228 229 230 !This section has been created automatically by the script Abilint (TD). 231 !Do not modify the following lines by hand. 232 #undef ABI_FUNC 233 #define ABI_FUNC 'pawcprj_free' 234 !End of the abilint section 235 236 implicit none 237 238 !Arguments ------------------------------------ 239 !scalars 240 !arrays 241 type(pawcprj_type),intent(inout) :: cprj(:,:) 242 243 !Local variables------------------------------- 244 !scalars 245 integer :: ii,jj,n1dim,n2dim 246 247 ! ************************************************************************* 248 249 n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2) 250 251 do jj=1,n2dim 252 do ii=1,n1dim 253 if (allocated(cprj(ii,jj)%cp)) then 254 LIBPAW_DEALLOCATE(cprj(ii,jj)%cp) 255 end if 256 if (allocated(cprj(ii,jj)%dcp)) then 257 LIBPAW_DEALLOCATE(cprj(ii,jj)%dcp) 258 end if 259 end do 260 end do 261 262 end subroutine pawcprj_free
m_pawcprj/pawcprj_gather_spin [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_gather_spin
FUNCTION
INPUTS
cprj(:,:)=the input cprj datastructure n2size=number of cprj datastructures to be gathered (second dim) nspinor : number of spinorial component (on current proc) nspinortot : total number of spinorial component
OUTPUT
cprj_gat(:,:) = the cprj containing all nspinor componants
NOTES
The cprj has been built like the following: loop on nsppol loop on k point loop over band or block of band These quantities were build only if treated by the current proc the inner quantities being nspinor
PARENTS
energy,pawmkrhoij
CHILDREN
xmpi_sum
SOURCE
2760 subroutine pawcprj_gather_spin(cprj,cprj_gat,natom,n2size,nspinor,nspinortot,& 2761 & spaceComm_spin,ierr) 2762 2763 2764 !This section has been created automatically by the script Abilint (TD). 2765 !Do not modify the following lines by hand. 2766 #undef ABI_FUNC 2767 #define ABI_FUNC 'pawcprj_gather_spin' 2768 !End of the abilint section 2769 2770 implicit none 2771 2772 !Arguments ------------------------------------ 2773 !scalars 2774 integer,intent(in) :: natom,nspinor,nspinortot,n2size 2775 integer,intent(in) :: spaceComm_spin 2776 integer,intent(out) :: ierr 2777 !arrays 2778 type(pawcprj_type),intent(in) :: cprj(:,:) 2779 type(pawcprj_type),intent(inout) :: cprj_gat(:,:) 2780 2781 !Local variables------------------------------- 2782 !scalars 2783 integer :: i1,iatom,ibsp,icpgr,ilmn,isp,ispinor,jj,lmndim,n2dim,n2dim_gat,ncpgr 2784 character(len=100) :: msg 2785 !arrays 2786 integer :: nlmn(natom) 2787 real(dp),allocatable :: buffer1(:),buffer2(:) 2788 2789 ! ************************************************************************* 2790 2791 n2dim =size(cprj,dim=2) 2792 n2dim_gat=size(cprj_gat,dim=2) 2793 if (n2dim_gat/=(nspinortot/nspinor)*n2dim) then 2794 msg='wrong dims (pawcprj_gather_spin)!' 2795 MSG_BUG(msg) 2796 end if 2797 2798 do iatom=1,natom 2799 nlmn(iatom)=size(cprj(iatom,1)%cp(1,:)) 2800 end do 2801 ncpgr=cprj(1,1)%ncpgr 2802 lmndim=2*n2size*sum(nlmn(1:natom))*(1+ncpgr) 2803 LIBPAW_ALLOCATE(buffer1,(lmndim)) 2804 LIBPAW_ALLOCATE(buffer2,(lmndim*nspinortot)) 2805 2806 isp=0;ibsp=0 2807 jj=1 2808 do i1=1,n2size 2809 isp=isp+1 2810 do iatom=1,natom 2811 do ilmn=1,nlmn(iatom) 2812 buffer1(jj:jj+1)=cprj(iatom,isp)%cp(1:2,ilmn) 2813 jj=jj+2 2814 end do 2815 if (ncpgr>0) then 2816 do ilmn=1,nlmn(iatom) 2817 do icpgr=1,ncpgr 2818 buffer1(jj:jj+1)=cprj(iatom,isp)%dcp(1:2,icpgr,ilmn) 2819 jj=jj+2 2820 end do 2821 end do 2822 end if 2823 end do 2824 end do 2825 2826 call xmpi_allgather(buffer1,lmndim,buffer2,spaceComm_spin,ierr) 2827 2828 jj=1 2829 do ispinor=1,nspinortot 2830 do i1 =1,n2size 2831 ibsp=(i1-1)*nspinortot + ispinor 2832 do iatom=1,natom 2833 do ilmn=1,nlmn(iatom) 2834 cprj_gat(iatom,ibsp)%cp(1:2,ilmn)=buffer2(jj:jj+1) 2835 jj=jj+2 2836 end do 2837 if (ncpgr>0) then 2838 do ilmn=1,nlmn(iatom) 2839 do icpgr=1,ncpgr 2840 cprj_gat(iatom,ibsp)%dcp(1:2,icpgr,ilmn)=buffer2(jj:jj+1) 2841 jj=jj+2 2842 end do 2843 end do 2844 end if 2845 end do 2846 end do 2847 end do 2848 2849 LIBPAW_DEALLOCATE(buffer1) 2850 LIBPAW_DEALLOCATE(buffer2) 2851 2852 end subroutine pawcprj_gather_spin
m_pawcprj/pawcprj_get [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_get
FUNCTION
Read the cprj_k for a given k-point from memory in cprj or from a temporary file
INPUTS
atind(natom)=index table for atoms (see iorder below) cprj(dimcp,nspinor*mband*mkmem*nsppol)=input cprj (used if mkmem/=0) dimcp=first dimension of cprj_k,cprj arrays (1 or natom) iband1=index of first band in cprj ibg=shift in cprj array to locate current k-point [icpgr]= (optional argument) if present, only component icpgr of input cprj gradient is copied into output cprj Not used if cprj(:,:)%ncpgr<icpgr (mkmem>0) or ncpgr(optional)<icpgr (mkmem=0) ikpt=index of current k-point (only needed for the parallel distribution) iorder=0 if cprj ordering does not change during reading 1 if cprj ordering changes during reading, depending on content of atind array: - if atind=atindx (type-sorted=>unsorted) - if atind=atindx1 (unsorted=>type-sorted) isppol=index of current spin component mband=maximum number of bands mkmem=number of k points which can fit in memory; set to 0 if use disk [mpi_comm]=(optional argument) MPI communicator over (k-pts,bands,spins) Must be used in association with proc_distrb argument natom=number of atoms in cell nband=number of bands to import (usually 1 or nband_k) nband_k=total number of bands for this k-point [ncpgr]=(optional argument) second dimension of cprj%dcp(2,ncpgr,nlmn stored in memory (mkmem>0) or present on disk (mkmem=0)) needed only when optional argument icpgr is present nspinor=number of spinorial components of the wavefunctions (on current proc) nsppol=1 for unpolarized, 2 for spin-polarized [proc_distrb(nkpt,nband,nsppol)]=(optional argument) processor distribution Describe how cprj datastructures are distributed over processors When present, mpicomm argument must be also present uncp=unit number for cprj data (used if mkmem=0)
OUTPUT
cprj_k(dimcp,nspinor*nband) <type(pawcprj_type)>= output cprj datastructure
PARENTS
berryphase_new,cgwf,chern_number,datafordmft,dfpt_nstpaw,dfpt_vtowfk dfpt_wfkfermi,dotprod_set_cgcprj,dotprodm_sumdiag_cgcprj,extrapwf fock2ACE,forstrnps,m_plowannier,make_grad_berry,optics_paw optics_paw_core,pawmkrhoij,posdoppler,rf2_init,smatrix_pawinit suscep_stat,wf_mixing
CHILDREN
xmpi_sum
SOURCE
1201 subroutine pawcprj_get(atind,cprj_k,cprj,dimcp,iband1,ibg,ikpt,iorder,isppol,mband,& 1202 & mkmem,natom,nband,nband_k,nspinor,nsppol,uncp,& 1203 & icpgr,ncpgr,mpicomm,proc_distrb) ! optionals arguments 1204 1205 1206 !This section has been created automatically by the script Abilint (TD). 1207 !Do not modify the following lines by hand. 1208 #undef ABI_FUNC 1209 #define ABI_FUNC 'pawcprj_get' 1210 !End of the abilint section 1211 1212 implicit none 1213 1214 !Arguments ------------------------------------ 1215 !scalars 1216 integer,intent(in) :: dimcp,iband1,ibg,ikpt,iorder,isppol,mband,mkmem,natom 1217 integer,intent(in) :: nband,nband_k,nspinor,nsppol,uncp 1218 integer,intent(in),optional :: icpgr,mpicomm,ncpgr 1219 !arrays 1220 integer,intent(in) :: atind(natom) 1221 integer,intent(in),optional :: proc_distrb(:,:,:) 1222 type(pawcprj_type),intent(in) :: cprj(dimcp,nspinor*mband*mkmem*nsppol) 1223 type(pawcprj_type),intent(inout) :: cprj_k(dimcp,nspinor*nband) 1224 1225 !Local variables------------------------------- 1226 !scalars 1227 integer :: iatm,iatom,ib,ibsp,icpgr_,isp,ispinor,jband,me,nband0,ncpgr_ 1228 logical :: has_distrb,has_icpgr 1229 character(len=500) :: msg 1230 !arrays 1231 real(dp),allocatable :: tmp(:,:,:) 1232 1233 ! ************************************************************************* 1234 1235 ncpgr_=cprj_k(1,1)%ncpgr;if (present(ncpgr)) ncpgr_=ncpgr 1236 icpgr_=-1;if(present(icpgr)) icpgr_=icpgr 1237 has_icpgr=(icpgr_>0.and.icpgr_<=ncpgr_) 1238 if (present(icpgr).and.(.not.present(ncpgr))) then 1239 msg='ncpgr must be present when icpgr is present (pawcprj_get)!' 1240 MSG_BUG(msg) 1241 end if 1242 if (has_icpgr.and.cprj_k(1,1)%ncpgr<1) then 1243 msg='cprj_k%ncpgr not consistent with icpgr (pawcprj_get)!' 1244 MSG_BUG(msg) 1245 end if 1246 1247 !MPI data 1248 has_distrb=present(proc_distrb) 1249 if (has_distrb) then 1250 if (.not.present(mpicomm)) then 1251 msg='mpicomm must be present when proc_distrb is present (pawcprj_get)!' 1252 MSG_BUG(msg) 1253 end if 1254 me=xmpi_comm_rank(mpicomm) 1255 end if 1256 1257 if (mkmem==0) then 1258 1259 if (iband1==1) then 1260 read(uncp) nband0 1261 if (nband_k/=nband0) then 1262 msg='_PAW file was not created with the right options (pawcprj_get)!' 1263 MSG_BUG(msg) 1264 end if 1265 end if 1266 1267 isp=0;jband=iband1-1 1268 do ib=1,nband 1269 jband=jband+1 1270 if (has_distrb) then 1271 if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then 1272 isp=isp+nspinor 1273 cycle 1274 end if 1275 end if 1276 do ispinor=1,nspinor 1277 isp=isp+1 1278 if (iorder==0) then 1279 if (ncpgr_==0) then 1280 do iatom=1,dimcp 1281 read(uncp) cprj_k(iatom,isp)%cp(:,:) 1282 end do 1283 else 1284 if (has_icpgr) then 1285 do iatom=1,dimcp 1286 LIBPAW_ALLOCATE(tmp,(2,ncpgr_,cprj_k(iatom,1)%nlmn)) 1287 read(uncp) cprj_k(iatom,isp)%cp(:,:),tmp(:,:,:) 1288 cprj_k(iatom,isp)%dcp(:,1,:)=tmp(:,icpgr_,:) 1289 LIBPAW_DEALLOCATE(tmp) 1290 end do 1291 else 1292 do iatom=1,dimcp 1293 read(uncp) cprj_k(iatom,isp)%cp(:,:),cprj_k(iatom,isp)%dcp(:,:,:) 1294 end do 1295 end if 1296 end if 1297 else 1298 if (ncpgr_==0) then 1299 do iatom=1,dimcp 1300 iatm=min(atind(iatom),dimcp) 1301 read(uncp) cprj_k(iatm,isp)%cp(:,:) 1302 end do 1303 else 1304 if (has_icpgr) then 1305 do iatom=1,dimcp 1306 iatm=min(atind(iatom),dimcp) 1307 LIBPAW_ALLOCATE(tmp,(2,ncpgr_,cprj_k(iatm,1)%nlmn)) 1308 read(uncp) cprj_k(iatm,isp)%cp(:,:),tmp(:,:,:) 1309 cprj_k(iatm,isp)%dcp(:,1,:)=tmp(:,icpgr_,:) 1310 LIBPAW_DEALLOCATE(tmp) 1311 end do 1312 else 1313 do iatom=1,dimcp 1314 iatm=min(atind(iatom),dimcp) 1315 read(uncp) cprj_k(iatm,isp)%cp(:,:),cprj_k(iatm,isp)%dcp(:,:,:) 1316 end do 1317 end if 1318 end if 1319 end if 1320 end do 1321 end do 1322 1323 else 1324 1325 isp=0;ibsp=ibg+nspinor*(iband1-1);jband=iband1-1 1326 do ib=1,nband 1327 jband=jband+1 1328 if (has_distrb) then 1329 if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then 1330 isp=isp+nspinor;ibsp=ibsp+nspinor 1331 cycle 1332 end if 1333 end if 1334 do ispinor=1,nspinor 1335 isp=isp+1;ibsp=ibsp+1 1336 if (iorder==0) then 1337 if (ncpgr_==0) then 1338 do iatom=1,dimcp 1339 cprj_k(iatom,isp)%cp(:,:)=cprj(iatom,ibsp)%cp(:,:) 1340 end do 1341 else 1342 if (has_icpgr) then 1343 do iatom=1,dimcp 1344 cprj_k(iatom,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:) 1345 cprj_k(iatom,isp)%dcp(:,1,:)=cprj(iatom,ibsp)%dcp(:,icpgr_,:) 1346 end do 1347 else 1348 do iatom=1,dimcp 1349 cprj_k(iatom,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:) 1350 cprj_k(iatom,isp)%dcp(:,:,:)=cprj(iatom,ibsp)%dcp(:,:,:) 1351 end do 1352 end if 1353 end if 1354 else 1355 if (ncpgr_==0) then 1356 do iatom=1,dimcp 1357 iatm=min(atind(iatom),dimcp) 1358 cprj_k(iatm,isp)%cp(:,:)=cprj(iatom,ibsp)%cp(:,:) 1359 end do 1360 else 1361 if (has_icpgr) then 1362 do iatom=1,dimcp 1363 iatm=min(atind(iatom),dimcp) 1364 cprj_k(iatm,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:) 1365 cprj_k(iatm,isp)%dcp(:,1,:)=cprj(iatom,ibsp)%dcp(:,icpgr_,:) 1366 end do 1367 else 1368 do iatom=1,dimcp 1369 iatm=min(atind(iatom),dimcp) 1370 cprj_k(iatm,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:) 1371 cprj_k(iatm,isp)%dcp(:,:,:)=cprj(iatom,ibsp)%dcp(:,:,:) 1372 end do 1373 end if 1374 end if 1375 end if 1376 end do 1377 end do 1378 1379 end if 1380 1381 end subroutine pawcprj_get
m_pawcprj/pawcprj_getdim [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_getdim
FUNCTION
Helper function returning the number of lmn components in the <p_{lmn}^i|\psi> for the i-th atom. Used to initialize the dimensioning array that is passed to the pawcprj_alloc routines when the pawcprj_type structure is allocated and initialized.
INPUTS
natom=number of atoms in the unit cell nattyp(ntypat)=number of atoms of each type ntypat=number of atom types typat(natom-= type of each atom Pawtab(ntypat)<pawtab_type>=PAW tabulated starting data. sort_mode(len=*)=String defining the sorting of the atoms in the Cprj arrays. Two modes are possible: -- "O[rdered]", if atoms are sorted by atom type. -- "R[andom]", if atoms are sorted randomly i.e. according the values of typat specified in the input file.
OUTPUT
dimcprj(natom)=Number of nlm elements in the <p_{lmn}^i|\psi> matrix elements for i=1,...,natom.
PARENTS
afterscfloop,berryphase_new,chern_number,dfpt_looppert,dfpt_scfcv extrapwf,forstr,getghc,initberry,m_fock,m_hamiltonian,mlwfovlp_qp outkss,scfcv,smatrix_pawinit,wf_mixing
CHILDREN
xmpi_sum
SOURCE
2890 subroutine pawcprj_getdim(dimcprj,natom,nattyp,ntypat,typat,Pawtab,sort_mode) 2891 2892 2893 !This section has been created automatically by the script Abilint (TD). 2894 !Do not modify the following lines by hand. 2895 #undef ABI_FUNC 2896 #define ABI_FUNC 'pawcprj_getdim' 2897 !End of the abilint section 2898 2899 implicit none 2900 2901 !Arguments ------------------------------------ 2902 integer,intent(in) :: natom,ntypat 2903 character(len=*),intent(in) :: sort_mode 2904 !arrays 2905 integer,intent(in) :: nattyp(:),typat(natom) 2906 integer,intent(inout) :: dimcprj(natom) 2907 type(Pawtab_type),intent(in) :: Pawtab(ntypat) 2908 2909 !Local variables------------------------------- 2910 integer :: iatom,itypat 2911 character(len=500) :: msg 2912 2913 ! ************************************************************************* 2914 2915 SELECT CASE (sort_mode(1:1)) 2916 2917 CASE ("o","O") ! Ordered by atom-type 2918 2919 iatom=0 2920 do itypat=1,ntypat 2921 dimcprj(iatom+1:iatom+nattyp(itypat))=Pawtab(itypat)%lmn_size 2922 iatom=iatom+nattyp(itypat) 2923 end do 2924 2925 CASE ("r","R") ! Randomly ordered (typat from input file) 2926 2927 do iatom=1,natom 2928 itypat=typat(iatom) 2929 dimcprj(iatom)=Pawtab(itypat)%lmn_size 2930 end do 2931 2932 CASE DEFAULT 2933 msg='Wrong value for sort_mode: '//TRIM(sort_mode) 2934 MSG_ERROR(msg) 2935 END SELECT 2936 2937 end subroutine pawcprj_getdim
m_pawcprj/pawcprj_lincom [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_lincom
FUNCTION
Compute a LINear COMbination of cprj datastructure: cprj_out(:,:) <--- Sum_i [ alpha_i . cprj_i(:,:) ] alpha_i are COMPLEX scalars
INPUTS
alpha(2,nn)= alpha COMPLEX factors cprj_in(:,:) <type(pawcprj_type)>= input cprj_in datastructure nn= number of cprj involved in the linear combination
OUTPUT
cprj_out(:,:) <type(pawcprj_type)>= output cprj_out datastructure
NOTES
cprj_in and cprj_out must be dimensionned as cprj_in(n1,n2*nn) and cprj_in(n1,n2)
PARENTS
extrapwf,getdc1,lincom_cgcprj,wf_mixing
CHILDREN
xmpi_sum
SOURCE
996 subroutine pawcprj_lincom(alpha,cprj_in,cprj_out,nn) 997 998 999 !This section has been created automatically by the script Abilint (TD). 1000 !Do not modify the following lines by hand. 1001 #undef ABI_FUNC 1002 #define ABI_FUNC 'pawcprj_lincom' 1003 !End of the abilint section 1004 1005 implicit none 1006 1007 !Arguments ------------------------------------ 1008 !scalars 1009 integer,intent(in) :: nn 1010 real(dp),intent(in) :: alpha(2,nn) 1011 !arrays 1012 type(pawcprj_type),intent(in) :: cprj_in(:,:) 1013 type(pawcprj_type),intent(inout) :: cprj_out(:,:) 1014 1015 !Local variables------------------------------- 1016 !scalars 1017 integer :: ii,in,jj,jn,kk,ll,n1in,n1out,n2in,n2out,ncpgrin,ncpgrout,nlmn 1018 real(dp) :: cp1,cp2 1019 character(len=500) :: msg 1020 1021 ! ************************************************************************* 1022 1023 n1in=size(cprj_in,dim=1);n1out=size(cprj_out,dim=1) 1024 n2in=size(cprj_in,dim=2);n2out=size(cprj_out,dim=2) 1025 ncpgrin=cprj_in(1,1)%ncpgr;ncpgrout=cprj_out(1,1)%ncpgr 1026 1027 msg = "" 1028 if (n1in/=n1out) msg = TRIM(msg)//"Bug in pawcprj_lincom: n1 wrong sizes!"//ch10 1029 if (n2in/=n2out*nn) msg = TRIM(msg)//"Bug in pawcprj_lincom: n2 wrong sizes!"//ch10 1030 if (ncpgrin/=ncpgrout) msg = TRIM(msg)//"Bug in pawcprj_lincom: ncpgr wrong sizes!"//ch10 1031 if (LEN_TRIM(msg) > 0) then 1032 MSG_ERROR(msg) 1033 end if 1034 1035 do jj=1,n2out 1036 do ii=1,n1out 1037 nlmn=cprj_in(ii,jj)%nlmn 1038 cprj_out(ii,jj)%nlmn=nlmn 1039 cprj_out(ii,jj)%cp(1:2,1:nlmn)=zero 1040 jn=jj 1041 do in=1,nn 1042 do kk=1,nlmn 1043 cp1=cprj_out(ii,jj)%cp(1,kk) & 1044 & +alpha(1,in)*cprj_in(ii,jn)%cp(1,kk)-alpha(2,in)*cprj_in(ii,jn)%cp(2,kk) 1045 cp2=cprj_out(ii,jj)%cp(2,kk) & 1046 & +alpha(1,in)*cprj_in(ii,jn)%cp(2,kk)+alpha(2,in)*cprj_in(ii,jn)%cp(1,kk) 1047 cprj_out(ii,jj)%cp(1,kk)=cp1 1048 cprj_out(ii,jj)%cp(2,kk)=cp2 1049 end do 1050 jn=jn+n2out 1051 end do 1052 end do 1053 end do 1054 1055 if (ncpgrin>0) then 1056 do jj=1,n2out 1057 do ii=1,n1out 1058 nlmn=cprj_in(ii,jj)%nlmn 1059 cprj_out(ii,jj)%dcp(1:2,1:ncpgrin,1:nlmn)=zero 1060 jn=jj 1061 do in=1,nn 1062 do kk=1,nlmn 1063 do ll=1,ncpgrin 1064 cp1=cprj_out(ii,jj)%dcp(1,ll,kk) & 1065 & +alpha(1,in)*cprj_in(ii,jn)%dcp(1,ll,kk) & 1066 & -alpha(2,in)*cprj_in(ii,jn)%dcp(2,ll,kk) 1067 cp2=cprj_out(ii,jj)%dcp(2,ll,kk) & 1068 & +alpha(1,in)*cprj_in(ii,jn)%dcp(2,ll,kk) & 1069 +alpha(2,in)*cprj_in(ii,jn)%dcp(1,ll,kk) 1070 cprj_out(ii,jj)%dcp(1,ll,kk)=cp1 1071 cprj_out(ii,jj)%dcp(2,ll,kk)=cp2 1072 end do 1073 end do 1074 jn=jn+n2out 1075 end do 1076 end do 1077 end do 1078 end if 1079 1080 end subroutine pawcprj_lincom
m_pawcprj/pawcprj_mpi_allgather [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_mpi_allgather
FUNCTION
Perform MPI_ALLGATHER on a pawcprj_type inside a MPI communicator.
INPUTS
cprj_loc= The cprj on the local proc being all-gathered natom=Number of atoms (size of first dimension of cprj_loc). n2dim=Size of the second dimension of cprj_loc. nlmn(natom)=Number of nlm partial waves for each atom. ncpgr = number of gradients in cprj_loc nproc=number of processors being gathered spaceComm=MPI Communicator. [rank_ordered]= optional, default=FALSE TRUE: second dimension of gathered datastructure is rank-ordered FALSE: second dimension of gathered datastructure is not rank-ordered
OUTPUT
cprj_gat=the gathered cprjs ierr=Error status.
PARENTS
berryphase_new,cgwf,optics_paw,optics_paw_core,suscep_stat
CHILDREN
xmpi_sum
SOURCE
2236 subroutine pawcprj_mpi_allgather(cprj_loc,cprj_gat,natom,n2dim,nlmn,ncpgr,nproc,spaceComm,ierr,& 2237 & rank_ordered) 2238 2239 2240 !This section has been created automatically by the script Abilint (TD). 2241 !Do not modify the following lines by hand. 2242 #undef ABI_FUNC 2243 #define ABI_FUNC 'pawcprj_mpi_allgather' 2244 !End of the abilint section 2245 2246 implicit none 2247 2248 !Arguments ------------------------------------ 2249 !scalars 2250 integer,intent(in) :: natom,n2dim,ncpgr,nproc,spaceComm 2251 integer,intent(out) :: ierr 2252 logical,optional,intent(in) :: rank_ordered 2253 !arrays 2254 integer,intent(in) :: nlmn(natom) 2255 type(pawcprj_type),intent(in) :: cprj_loc(:,:) 2256 type(pawcprj_type),intent(inout) :: cprj_gat(:,:) !vz_i 2257 2258 !Local variables------------------------------- 2259 !scalars 2260 integer :: iat,jj,t2dim,tcpgr,tg2dim,n1dim,nn 2261 integer :: ntotcp,ibuf,ipck,iproc 2262 logical :: rank_ordered_ 2263 character(len=100) :: msg 2264 !arrays 2265 real(dp),allocatable :: buffer_cpgr(:,:,:),buffer_cpgr_all(:,:,:) 2266 2267 ! ************************************************************************* 2268 2269 n1dim=0 2270 t2dim=0 2271 tg2dim=0 2272 tcpgr=0 2273 ierr=0 2274 2275 nn=size(nlmn,dim=1) 2276 n1dim=size(cprj_loc,dim=1) 2277 t2dim=size(cprj_loc,dim=2) 2278 tg2dim=size(cprj_gat,dim=2) 2279 tcpgr=cprj_loc(1,1)%ncpgr 2280 2281 if (nn/=n1dim) then 2282 msg='size mismatch in natom (pawcprj_mpi_allgather)!' 2283 MSG_BUG(msg) 2284 end if 2285 if (t2dim/=n2dim) then 2286 msg='size mismatch in dim=2 (pawcprj_mpi_allgather)!' 2287 MSG_BUG(msg) 2288 end if 2289 if (tg2dim/=n2dim*nproc) then 2290 msg='size mismatch in dim=2 (pawcprj_mpi_allgather)!' 2291 MSG_BUG(msg) 2292 end if 2293 if (tcpgr/=ncpgr) then 2294 msg='size mismatch in ncpgr (pawcprj_mpi_allgather)!' 2295 MSG_BUG(msg) 2296 end if 2297 2298 rank_ordered_=.false.;if(present(rank_ordered)) rank_ordered_=rank_ordered 2299 2300 ntotcp=n2dim*SUM(nlmn(:)) 2301 LIBPAW_ALLOCATE(buffer_cpgr,(2,1+ncpgr,ntotcp)) 2302 LIBPAW_ALLOCATE(buffer_cpgr_all,(2,1+ncpgr,nproc*ntotcp)) 2303 2304 !=== Pack cprj_loc ==== 2305 ipck=0 2306 do jj=1,n2dim 2307 do iat=1,natom 2308 nn=nlmn(iat) 2309 buffer_cpgr(:,1,ipck+1:ipck+nn)=cprj_loc(iat,jj)%cp(:,1:nn) 2310 if (ncpgr/=0) buffer_cpgr(:,2:1+ncpgr,ipck+1:ipck+nn)=cprj_loc(iat,jj)%dcp(:,:,1:nn) 2311 ipck=ipck+nn 2312 end do 2313 end do 2314 2315 !=== allgather data === 2316 call xmpi_allgather(buffer_cpgr,2*(ncpgr+1)*ntotcp,buffer_cpgr_all,spaceComm,ierr) 2317 2318 !=== unpack gathered data into cprj(natom,n2dim*nproc) 2319 !=== second dimension is rank-ordered if rank_ordered_=true 2320 ipck=0 2321 do iproc=1,nproc 2322 do jj=1,n2dim 2323 if (rank_ordered_) then 2324 ibuf=jj+(iproc-1)*n2dim 2325 else 2326 ibuf=iproc+(jj-1)*nproc 2327 end if 2328 do iat=1,natom 2329 nn=nlmn(iat) 2330 cprj_gat(iat,ibuf)%cp(:,1:nn)=buffer_cpgr_all(:,1,ipck+1:ipck+nn) 2331 if (ncpgr/=0) cprj_gat(iat,ibuf)%dcp(:,1:ncpgr,1:nn)=& 2332 & buffer_cpgr_all(:,2:1+ncpgr,ipck+1:ipck+nn) 2333 ipck=ipck+nn 2334 end do 2335 end do 2336 end do 2337 2338 !=== Clean up === 2339 LIBPAW_DEALLOCATE(buffer_cpgr) 2340 LIBPAW_DEALLOCATE(buffer_cpgr_all) 2341 2342 end subroutine pawcprj_mpi_allgather
m_pawcprj/pawcprj_mpi_exch [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_mpi_exch
FUNCTION
Exchange a pawcprj_type between two processors inside a MPI communicator.
INPUTS
natom=Number of atoms (size of first dimension of Cprj_send and Cprj_recv). n2dim=Size of the second dimension. nlmn(natom)=Number of nlm partial waves for each atom. Cprj_send= The datatype to be transmitted. receiver=ID of the receiver in spaceComm. sender=ID of the sender in spaceComm. spaceComm=MPI Communicator.
OUTPUT
ierr=Error status. Cprj_recv=The datatype copied on proc. receiver.
NOTES
If sender==receiver, Cprj_send is copied into Cprj_recv. It should be easy to avoid this additional copy in the calling routine.
PARENTS
outkss
CHILDREN
xmpi_sum
SOURCE
1755 subroutine pawcprj_mpi_exch(natom,n2dim,nlmn,ncpgr,Cprj_send,Cprj_recv,sender,receiver,spaceComm,ierr) 1756 1757 1758 !This section has been created automatically by the script Abilint (TD). 1759 !Do not modify the following lines by hand. 1760 #undef ABI_FUNC 1761 #define ABI_FUNC 'pawcprj_mpi_exch' 1762 !End of the abilint section 1763 1764 implicit none 1765 1766 !Arguments ------------------------------------ 1767 !scalars 1768 integer,intent(in) :: natom,n2dim,ncpgr 1769 integer,intent(in) :: sender,receiver,spaceComm 1770 integer,intent(out) :: ierr 1771 !arrays 1772 integer,intent(in) :: nlmn(natom) 1773 type(pawcprj_type),intent(in) :: Cprj_send(:,:) 1774 type(pawcprj_type),intent(inout) :: Cprj_recv(:,:) 1775 1776 !Local variables------------------------------- 1777 !scalars 1778 integer :: iat,jj,t2dim,tcpgr,n1dim,nn 1779 integer :: ntotcp,ipck,rank 1780 character(len=500) :: msg 1781 !arrays 1782 real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:) 1783 1784 ! ************************************************************************* 1785 1786 n1dim=0 1787 t2dim=0 1788 tcpgr=0 1789 ierr=0 1790 if (sender==receiver) then 1791 call pawcprj_copy(Cprj_send,Cprj_recv) 1792 return 1793 end if 1794 1795 rank = xmpi_comm_rank(spaceComm) 1796 1797 nn=size(nlmn,dim=1) 1798 if (rank==sender) then 1799 n1dim=size(Cprj_send,dim=1) 1800 t2dim=size(Cprj_send,dim=2) 1801 tcpgr=Cprj_send(1,1)%ncpgr 1802 end if 1803 if (rank==receiver) then 1804 n1dim=size(Cprj_recv,dim=1) 1805 t2dim=size(Cprj_recv,dim=2) 1806 tcpgr=Cprj_recv(1,1)%ncpgr 1807 end if 1808 if (rank/=sender.and.rank/=receiver) then 1809 write(msg,'(a,3i0)') & 1810 & 'rank is not equal to sender or receiver (pawcprj_mpi_exch): ',rank, sender, receiver 1811 MSG_BUG(msg) 1812 end if 1813 1814 ntotcp=n2dim*SUM(nlmn(:)) 1815 1816 LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp)) 1817 if (ncpgr/=0) then 1818 LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp)) 1819 end if 1820 1821 !=== Pack Cprj_send === 1822 if (rank==sender) then 1823 ipck=0 1824 do jj=1,n2dim 1825 do iat=1,natom 1826 nn=nlmn(iat) 1827 buffer_cp(:,ipck+1:ipck+nn)=Cprj_send(iat,jj)%cp(:,1:nn) 1828 if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=Cprj_send(iat,jj)%dcp(:,:,1:nn) 1829 ipck=ipck+nn 1830 end do 1831 end do 1832 end if 1833 1834 !=== Transmit data === 1835 call xmpi_exch(buffer_cp,2*ntotcp,sender,buffer_cp,receiver,spaceComm,ierr) 1836 if (ncpgr/=0) then 1837 call xmpi_exch(buffer_cpgr,2*ncpgr*ntotcp,sender,buffer_cpgr,receiver,spaceComm,ierr) 1838 end if 1839 1840 !=== UnPack buffers into Cprj_recv === 1841 if (rank==receiver) then 1842 ipck=0 1843 do jj=1,n2dim 1844 do iat=1,natom 1845 nn=nlmn(iat) 1846 Cprj_recv(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn) 1847 if (ncpgr/=0) Cprj_recv(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn) 1848 ipck=ipck+nn 1849 end do 1850 end do 1851 end if 1852 1853 LIBPAW_DEALLOCATE(buffer_cp) 1854 if (ncpgr/=0) then 1855 LIBPAW_DEALLOCATE(buffer_cpgr) 1856 end if 1857 1858 end subroutine pawcprj_mpi_exch
m_pawcprj/pawcprj_mpi_recv [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_mpi_recv
FUNCTION
Receive a pawcprj_type inside a MPI communicator.
INPUTS
natom=Number of atoms (size of first dimension of Cprj_in). n2dim=Size of the second dimension. nlmn(natom)=Number of nlm partial waves for each atom. ncpgr = number of gradients in cprj_in sender=ID of the sender in spaceComm. spaceComm=MPI Communicator.
OUTPUT
ierr=Error status. cprj_in=The datatype copied on proc. receiver.
NOTES
Perhaps in general it is more efficient to use pawcprj_mpi_exch but it is convenient for coding to have separate send and receive routines.
PARENTS
berryphase_new,posdoppler
CHILDREN
xmpi_sum
SOURCE
2015 subroutine pawcprj_mpi_recv(natom,n2dim,nlmn,ncpgr,cprj_in,sender,spaceComm,ierr) 2016 2017 2018 !This section has been created automatically by the script Abilint (TD). 2019 !Do not modify the following lines by hand. 2020 #undef ABI_FUNC 2021 #define ABI_FUNC 'pawcprj_mpi_recv' 2022 !End of the abilint section 2023 2024 implicit none 2025 2026 !Arguments ------------------------------------ 2027 !scalars 2028 integer,intent(in) :: natom,n2dim,ncpgr 2029 integer,intent(in) :: sender,spaceComm 2030 integer,intent(out) :: ierr 2031 !arrays 2032 integer,intent(in) :: nlmn(natom) 2033 type(pawcprj_type),intent(inout) :: cprj_in(:,:) 2034 2035 !Local variables------------------------------- 2036 !scalars 2037 integer :: iat,jj,t2dim,tcpgr,n1dim,nn 2038 integer :: ntotcp,ipck,tag 2039 character(len=100) :: msg 2040 !arrays 2041 real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:) 2042 2043 ! ************************************************************************* 2044 2045 n1dim=0 2046 t2dim=0 2047 tcpgr=0 2048 ierr=0 2049 2050 nn=size(nlmn,dim=1) 2051 n1dim=size(cprj_in,dim=1) 2052 t2dim=size(cprj_in,dim=2) 2053 tcpgr=cprj_in(1,1)%ncpgr 2054 2055 if (nn/=n1dim) then 2056 msg='size mismatch in natom (pawcprj_mpi_recv)!' 2057 MSG_BUG(msg) 2058 end if 2059 if (t2dim/=n2dim) then 2060 msg='size mismatch in dim=2 (pawcprj_mpi_recv)!' 2061 MSG_BUG(msg) 2062 end if 2063 if (tcpgr/=ncpgr) then 2064 msg='size mismatch in ncpgr (pawcprj_mpi_recv)!' 2065 MSG_BUG(msg) 2066 end if 2067 2068 ntotcp=n2dim*SUM(nlmn(:)) 2069 2070 LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp)) 2071 if (ncpgr/=0) then 2072 LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp)) 2073 end if 2074 2075 !=== Receive data === 2076 tag = 2*ntotcp 2077 call xmpi_recv(buffer_cp,sender,tag,spaceComm,ierr) 2078 if (ncpgr/=0) then 2079 tag=tag*ncpgr 2080 call xmpi_recv(buffer_cpgr,sender,tag,spaceComm,ierr) 2081 end if 2082 2083 !=== UnPack buffers into cprj_in === 2084 ipck=0 2085 do jj=1,n2dim 2086 do iat=1,natom 2087 nn=nlmn(iat) 2088 cprj_in(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn) 2089 if (ncpgr/=0) cprj_in(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn) 2090 ipck=ipck+nn 2091 end do 2092 end do 2093 2094 !=== Clean up === 2095 LIBPAW_DEALLOCATE(buffer_cp) 2096 if (ncpgr/=0) then 2097 LIBPAW_DEALLOCATE(buffer_cpgr) 2098 end if 2099 2100 end subroutine pawcprj_mpi_recv
m_pawcprj/pawcprj_mpi_send [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_mpi_send
FUNCTION
Send a pawcprj_type inside a MPI communicator.
INPUTS
natom=Number of atoms (size of first dimension of cprj_out). n2dim=Size of the second dimension. nlmn(natom)=Number of nlm partial waves for each atom. ncpgr = number of gradients in cprj_out cprj_out= The datatype to be transmitted. receiver=ID of the receiver in spaceComm. spaceComm=MPI Communicator.
OUTPUT
ierr=Error status.
NOTES
perhaps in general it is more efficient to use pawcprj_mpi_exch but it is convenient for coding to have separate send and recieve routines.
PARENTS
berryphase_new,posdoppler
CHILDREN
xmpi_sum
SOURCE
1894 subroutine pawcprj_mpi_send(natom,n2dim,nlmn,ncpgr,cprj_out,receiver,spaceComm,ierr) 1895 1896 1897 !This section has been created automatically by the script Abilint (TD). 1898 !Do not modify the following lines by hand. 1899 #undef ABI_FUNC 1900 #define ABI_FUNC 'pawcprj_mpi_send' 1901 !End of the abilint section 1902 1903 implicit none 1904 1905 !Arguments ------------------------------------ 1906 !scalars 1907 integer,intent(in) :: natom,n2dim,ncpgr 1908 integer,intent(in) :: receiver,spaceComm 1909 integer,intent(out) :: ierr 1910 !arrays 1911 integer,intent(in) :: nlmn(natom) 1912 type(pawcprj_type),intent(in) :: cprj_out(:,:) 1913 1914 !Local variables------------------------------- 1915 !scalars 1916 integer :: iat,jj,t2dim,tcpgr,n1dim,nn 1917 integer :: ntotcp,ipck,tag 1918 character(len=100) :: msg 1919 !arrays 1920 real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:) 1921 1922 ! ************************************************************************* 1923 1924 n1dim=0 1925 t2dim=0 1926 tcpgr=0 1927 ierr=0 1928 1929 nn=size(nlmn,dim=1) 1930 n1dim=size(cprj_out,dim=1) 1931 t2dim=size(cprj_out,dim=2) 1932 tcpgr=cprj_out(1,1)%ncpgr 1933 1934 if (nn/=n1dim) then 1935 msg='size mismatch in natom (pawcprj_mpi_send)!' 1936 MSG_BUG(msg) 1937 end if 1938 if (t2dim/=n2dim) then 1939 msg='size mismatch in dim=2 (pawcprj_mpi_send)!' 1940 MSG_BUG(msg) 1941 end if 1942 if (tcpgr/=ncpgr) then 1943 msg='size mismatch in ncpgr (pawcprj_mpi_send)!' 1944 MSG_BUG(msg) 1945 end if 1946 1947 ntotcp=n2dim*SUM(nlmn(:)) 1948 1949 LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp)) 1950 if (ncpgr/=0) then 1951 LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp)) 1952 end if 1953 1954 !=== Pack cprj_out ==== 1955 ipck=0 1956 do jj=1,n2dim 1957 do iat=1,natom 1958 nn=nlmn(iat) 1959 buffer_cp(:,ipck+1:ipck+nn)=cprj_out(iat,jj)%cp(:,1:nn) 1960 if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=cprj_out(iat,jj)%dcp(:,:,1:nn) 1961 ipck=ipck+nn 1962 end do 1963 end do 1964 1965 !=== Transmit data === 1966 tag = 2*ntotcp 1967 call xmpi_send(buffer_cp,receiver,tag,spaceComm,ierr) 1968 if (ncpgr/=0) then 1969 tag=tag*ncpgr 1970 call xmpi_send(buffer_cpgr,receiver,tag,spaceComm,ierr) 1971 end if 1972 1973 !=== Clean up === 1974 LIBPAW_DEALLOCATE(buffer_cp) 1975 if (ncpgr/=0) then 1976 LIBPAW_DEALLOCATE(buffer_cpgr) 1977 end if 1978 1979 end subroutine pawcprj_mpi_send
m_pawcprj/pawcprj_mpi_sum [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_mpi_sum
FUNCTION
Perform MPI_SUM on a pawcprj_type inside a MPI communicator.
INPUTS
spaceComm=MPI Communicator.
SIDE EFFECTS
cprj=the cprj datastructure ierr=Error status.
PARENTS
ctocprj
CHILDREN
xmpi_sum
SOURCE
2127 subroutine pawcprj_mpi_sum(cprj,spaceComm,ierr) 2128 2129 2130 !This section has been created automatically by the script Abilint (TD). 2131 !Do not modify the following lines by hand. 2132 #undef ABI_FUNC 2133 #define ABI_FUNC 'pawcprj_mpi_sum' 2134 !End of the abilint section 2135 2136 implicit none 2137 2138 !Arguments ------------------------------------ 2139 !scalars 2140 integer,intent(in) :: spaceComm 2141 integer,intent(out) :: ierr 2142 !arrays 2143 type(pawcprj_type),intent(inout) :: cprj(:,:) 2144 2145 !Local variables------------------------------- 2146 !scalars 2147 integer,parameter :: maxBytes=100*1024*1024 ! 100 MBytes 2148 integer :: ii,ipck,jj,ncpgr,nlmn,nn,n1dim,n2dim,n2dim1,n2dim2,sizeBytes,step 2149 logical,parameter :: save_memory=.true. 2150 !arrays 2151 real(dp),allocatable :: buffer_cprj(:,:,:) 2152 2153 ! ************************************************************************* 2154 2155 if (xmpi_comm_size(spaceComm)<2) return 2156 2157 n1dim=size(cprj,1);n2dim=size(cprj,2) 2158 nlmn=sum(cprj(:,:)%nlmn) 2159 ncpgr=maxval(cprj(:,:)%ncpgr) 2160 2161 step=n2dim 2162 if (save_memory) then 2163 sizeBytes=2*(1+ncpgr)*nlmn *8 2164 step=n2dim/max(1,sizeBytes/maxBytes) 2165 if (step==0) step=1 2166 end if 2167 2168 do n2dim1=1,n2dim,step 2169 2170 n2dim2=min(n2dim1+step-1,n2dim) 2171 nlmn=sum(cprj(:,n2dim1:n2dim2)%nlmn) 2172 LIBPAW_ALLOCATE(buffer_cprj,(2,1+ncpgr,nlmn)) 2173 2174 ipck=0 ; buffer_cprj=zero 2175 do jj=n2dim1,n2dim2 2176 do ii=1,n1dim 2177 nn=cprj(ii,jj)%nlmn 2178 buffer_cprj(:,1,ipck+1:ipck+nn)=cprj(ii,jj)%cp(:,1:nn) 2179 if (cprj(ii,jj)%ncpgr/=0) buffer_cprj(:,2:1+ncpgr,ipck+1:ipck+nn)=cprj(ii,jj)%dcp(:,1:ncpgr,1:nn) 2180 ipck=ipck+nn 2181 end do 2182 end do 2183 2184 call xmpi_sum(buffer_cprj,spaceComm,ierr) 2185 2186 ipck=0 2187 do jj=n2dim1,n2dim2 2188 do ii=1,n1dim 2189 nn=cprj(ii,jj)%nlmn 2190 cprj(ii,jj)%cp(:,1:nn)=buffer_cprj(:,1,ipck+1:ipck+nn) 2191 if (cprj(ii,jj)%ncpgr/=0) cprj(ii,jj)%dcp(:,1:ncpgr,1:nn)=buffer_cprj(:,2:1+ncpgr,ipck+1:ipck+nn) 2192 ipck=ipck+nn 2193 end do 2194 end do 2195 2196 LIBPAW_DEALLOCATE(buffer_cprj) 2197 2198 end do 2199 2200 end subroutine pawcprj_mpi_sum
m_pawcprj/pawcprj_output [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_output
FUNCTION
Output a cprj. Useful for debugging.
INPUTS
cprj(:,:) <type(pawcprj_type)>= cprj datastructure
OUTPUT
PARENTS
CHILDREN
xmpi_sum
SOURCE
1104 subroutine pawcprj_output(cprj) 1105 1106 1107 !This section has been created automatically by the script Abilint (TD). 1108 !Do not modify the following lines by hand. 1109 #undef ABI_FUNC 1110 #define ABI_FUNC 'pawcprj_output' 1111 !End of the abilint section 1112 1113 implicit none 1114 1115 !Arguments ------------------------------------ 1116 !scalar 1117 !arrays 1118 type(pawcprj_type),intent(in) :: cprj(:,:) 1119 1120 !Local variables------------------------------- 1121 !scalar 1122 integer :: ii,jj,kk,nlmn,n1dim,n2dim 1123 1124 ! ************************************************************************* 1125 1126 n1dim=size(cprj,dim=1) 1127 n2dim=size(cprj,dim=2) 1128 1129 write(std_out,'(a)')' pawcprj_output ' 1130 1131 do jj=1,n2dim 1132 do ii=1,n1dim 1133 write(std_out,'(a,i4,a,i4)')'atom ',ii,' band*k ',jj 1134 nlmn=cprj(ii,jj)%nlmn 1135 do kk=1,nlmn 1136 write(std_out,'(2f12.8)')cprj(ii,jj)%cp(1,kk),cprj(ii,jj)%cp(2,kk) 1137 end do 1138 end do 1139 end do 1140 1141 end subroutine pawcprj_output
m_pawcprj/pawcprj_put [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_put
FUNCTION
Write cprj_k for a given set of (n,k) into memory in cprj, or into a temporary file
INPUTS
atind(natom)=index table for atoms (see iorder below) cprj_k(dimcp,nspinor*nband) <type(pawcprj_type)>= input cprj datastructure dimcp=first dimension of cprj_k,cprjnk arrays (1 or natom) iband1=index of first band in cprj ibg=shift in cprj array to locate current k-point ikpt=index of current k-point (only needed for the parallel distribution) iorder=0 if cprj ordering does not change during reading 1 if cprj ordering changes during writing, depending on content of atind array: - if atind=atindx (type-sorted->unsorted) - if atind=atindx1 (unsorted->type-sorted) isppol=index of current spin component mband=maximum number of bands mkmem=number of k points which can fit in memory; set to 0 if use disk [mpi_comm]=(optional argument) MPI communicator over (k-pts,bands,spins) Must be used in association with proc_distrb argument [mpi_comm_band]=(optional argument) MPI communicator over bands Must be used in association with proc_distrb argument natom=number of atoms in cell nband=number of bands to export (usually 1, nband_k or nblockbd) nband_k=total number of bands for this k-point nlmn(dimcp)=array of dimensions of cprj_k,cprjnk datastructures nspinor=number of spinorial components of the wavefunctions (on current proc) nsppol=1 for unpolarized, 2 for spin-polarized [proc_distrb(nkpt,nband,nsppol)]=(optional argument) processor distribution Describe how cprj datastructures are distributed over processors When present, mpicomm argument must be also present [to_be_gathered]=(optional argument) TRUE if cprj_k arrays have to be gathered between procs (band-fft parallelism only) uncp=unit number for cprj data (used if mkmem=0)
SIDE EFFECTS
cprj(dimcp,nspinor*mband*mkmem*nsppol)=output cprj (used if mkmem/=0)
PARENTS
berryphase_new,cgwf,ctocprj,dfpt_vtowfk,extrapwf,vtowfk,wf_mixing
CHILDREN
xmpi_sum
SOURCE
1435 subroutine pawcprj_put(atind,cprj_k,cprj,dimcp,iband1,ibg,ikpt,iorder,isppol,mband,& 1436 & mkmem,natom,nband,nband_k,nlmn,nspinor,nsppol,uncp,& 1437 & mpicomm,mpi_comm_band,proc_distrb,to_be_gathered) ! Optional arguments 1438 1439 1440 !This section has been created automatically by the script Abilint (TD). 1441 !Do not modify the following lines by hand. 1442 #undef ABI_FUNC 1443 #define ABI_FUNC 'pawcprj_put' 1444 !End of the abilint section 1445 1446 implicit none 1447 1448 !Arguments ------------------------------------ 1449 !scalars 1450 integer,intent(in) :: iband1,ibg,ikpt,iorder,isppol,dimcp,mband,mkmem 1451 integer,intent(in) :: natom,nband,nband_k,nspinor,nsppol,uncp 1452 integer,intent(in),optional :: mpicomm,mpi_comm_band 1453 logical,optional,intent(in) :: to_be_gathered 1454 !arrays 1455 integer,intent(in) :: atind(natom),nlmn(dimcp) 1456 integer,intent(in),optional :: proc_distrb(:,:,:) 1457 type(pawcprj_type),intent(inout) :: cprj(dimcp,nspinor*mband*mkmem*nsppol) 1458 type(pawcprj_type),intent(in) :: cprj_k(dimcp,nspinor*nband) 1459 1460 !Local variables------------------------------- 1461 !scalars 1462 integer :: iatm,iatom,iband,ibsp,icpgr,ierr,ii,ilmn,isp,ispinor,jband,jj 1463 integer :: lmndim,me,ncpgr,nproc_band 1464 logical :: has_distrb,to_be_gathered_ 1465 character(len=500) :: msg 1466 !arrays 1467 real(dp),allocatable :: buffer1(:),buffer2(:) 1468 ! ************************************************************************* 1469 1470 ncpgr=cprj_k(1,1)%ncpgr 1471 to_be_gathered_=.false.;if (present(to_be_gathered)) to_be_gathered_=to_be_gathered 1472 1473 !MPI data 1474 nproc_band=1;if (present(mpi_comm_band)) nproc_band=xmpi_comm_size(mpi_comm_band) 1475 has_distrb=present(proc_distrb) 1476 if (has_distrb) then 1477 if (.not.present(mpicomm)) then 1478 msg='mpicomm must be present when proc_distrb is present (pawcprj_put)!' 1479 MSG_BUG(msg) 1480 end if 1481 me=xmpi_comm_rank(mpicomm) 1482 end if 1483 1484 if (nproc_band==1.or.(.not.to_be_gathered_)) then 1485 1486 if (mkmem==0) then 1487 1488 if (iband1==1) write(uncp) nband_k 1489 1490 isp=0;jband=iband1-1 1491 do iband=1,nband 1492 jband=jband+1 1493 if (has_distrb) then 1494 if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then 1495 isp=isp+nspinor 1496 cycle 1497 end if 1498 end if 1499 do ispinor=1,nspinor 1500 isp=isp+1 1501 if (iorder==0) then 1502 do iatom=1,dimcp 1503 if (ncpgr==0) then 1504 write(uncp) cprj_k(iatom,isp)%cp(:,:) 1505 else 1506 write(uncp) cprj_k(iatom,isp)%cp(:,:),cprj_k(iatom,isp)%dcp(:,:,:) 1507 end if 1508 end do 1509 else 1510 do iatom=1,dimcp 1511 iatm=min(atind(iatom),dimcp) 1512 if (ncpgr==0) then 1513 write(uncp) cprj_k(iatm,isp)%cp(:,:) 1514 else 1515 write(uncp) cprj_k(iatm,isp)%cp(:,:),cprj_k(iatm,isp)%dcp(:,:,:) 1516 end if 1517 end do 1518 end if 1519 end do 1520 end do 1521 1522 else 1523 1524 isp=0;ibsp=ibg+nspinor*(iband1-1);jband=iband1-1 1525 do iband=1,nband 1526 jband=jband+1 1527 if (has_distrb) then 1528 if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then 1529 isp=isp+nspinor;ibsp=ibsp+nspinor 1530 cycle 1531 end if 1532 end if 1533 do ispinor=1,nspinor 1534 isp=isp+1;ibsp=ibsp+1 1535 if (iorder==0) then 1536 do iatom=1,dimcp 1537 cprj(iatom,ibsp)%cp(:,:)=cprj_k(iatom,isp)%cp(:,:) 1538 if (ncpgr>0) cprj(iatom,ibsp)%dcp(:,:,:)=cprj_k(iatom,isp)%dcp(:,:,:) 1539 end do 1540 else 1541 do iatom=1,dimcp 1542 iatm=min(atind(iatom),dimcp) 1543 cprj(iatom,ibsp)%cp(:,:)=cprj_k(iatm,isp)%cp(:,:) 1544 if (ncpgr>0) cprj(iatom,ibsp)%dcp(:,:,:)=cprj_k(iatm,isp)%dcp(:,:,:) 1545 end do 1546 end if 1547 end do 1548 end do 1549 1550 end if 1551 1552 else ! np_band>1 1553 1554 lmndim=2*sum(nlmn(1:dimcp))*(1+ncpgr)*nspinor 1555 LIBPAW_ALLOCATE(buffer1,(lmndim)) 1556 LIBPAW_ALLOCATE(buffer2,(lmndim*nproc_band)) 1557 isp=0;ibsp=ibg+nspinor*(iband1-1) 1558 do iband=1,nband ! must be nblockbd for band-fft parallelism 1559 jj=1 1560 do ispinor=1,nspinor 1561 isp=isp+1 1562 do iatom=1,dimcp 1563 if (iorder==0) then 1564 iatm=iatom 1565 else 1566 iatm=min(atind(iatom),dimcp) 1567 end if 1568 do ilmn=1,nlmn(iatm) 1569 buffer1(jj:jj+1)=cprj_k(iatm,isp)%cp(1:2,ilmn) 1570 jj=jj+2 1571 end do 1572 if (ncpgr>0) then 1573 do ilmn=1,nlmn(iatm) 1574 do icpgr=1,ncpgr 1575 buffer1(jj:jj+1)=cprj_k(iatm,isp)%dcp(1:2,icpgr,ilmn) 1576 jj=jj+2 1577 end do 1578 end do 1579 end if 1580 end do !iatom 1581 end do !ispinor 1582 call xmpi_allgather(buffer1,lmndim,buffer2,mpi_comm_band,ierr) 1583 jj=1 1584 do ii=1,nproc_band 1585 do ispinor=1,nspinor 1586 ibsp=ibsp+1 1587 do iatom=1,dimcp 1588 if (iorder==0) then 1589 iatm=iatom 1590 else 1591 iatm=min(atind(iatom),dimcp) 1592 end if 1593 do ilmn=1,nlmn(iatm) 1594 cprj(iatom,ibsp)%cp(1:2,ilmn)=buffer2(jj:jj+1) 1595 jj=jj+2 1596 end do 1597 if (ncpgr>0) then 1598 do ilmn=1,nlmn(iatm) 1599 do icpgr=1,ncpgr 1600 cprj(iatom,ibsp)%dcp(1:2,icpgr,ilmn)=buffer2(jj:jj+1) 1601 jj=jj+2 1602 end do 1603 end do 1604 end if 1605 end do !iatom 1606 end do !ispinor 1607 end do !ii=1,nproc_band 1608 end do !iband 1609 LIBPAW_DEALLOCATE(buffer1) 1610 LIBPAW_DEALLOCATE(buffer2) 1611 1612 end if ! mode_para=b, nband 1613 1614 end subroutine pawcprj_put
m_pawcprj/pawcprj_reorder [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_reorder
FUNCTION
Change the order of a cprj datastructure From unsorted cprj to atom-sorted cprj (atm_indx=atindx) From atom-sorted cprj to unsorted cprj (atm_indx=atindx1)
INPUTS
atm_indx(natom)=index table for atoms From unsorted cprj to atom-sorted cprj (atm_indx=atindx) From atom-sorted cprj to unsorted cprj (atm_indx=atindx1)
OUTPUT
SIDE EFFECTS
cprj(:,:) <type(pawcprj_type)>= cprj datastructure
PARENTS
fock2ACE,forstrnps,ks_ddiago,scfcv
CHILDREN
xmpi_sum
SOURCE
1646 subroutine pawcprj_reorder(cprj,atm_indx) 1647 1648 1649 !This section has been created automatically by the script Abilint (TD). 1650 !Do not modify the following lines by hand. 1651 #undef ABI_FUNC 1652 #define ABI_FUNC 'pawcprj_reorder' 1653 !End of the abilint section 1654 1655 implicit none 1656 1657 !Arguments ------------------------------------ 1658 !scalars 1659 !arrays 1660 integer,intent(in) :: atm_indx(:) 1661 type(pawcprj_type),intent(inout) :: cprj(:,:) 1662 1663 !Local variables------------------------------- 1664 !scalars 1665 integer :: iexit,ii,jj,kk,n1atindx,n1cprj,n2cprj,ncpgr 1666 character(len=100) :: msg 1667 !arrays 1668 integer,allocatable :: nlmn(:) 1669 type(pawcprj_type),allocatable :: cprj_tmp(:,:) 1670 1671 ! ************************************************************************* 1672 1673 n1cprj=size(cprj,dim=1);n2cprj=size(cprj,dim=2) 1674 n1atindx=size(atm_indx,dim=1) 1675 if (n1cprj==0.or.n2cprj==0.or.n1atindx<=1) return 1676 1677 if (n1cprj/=n1atindx) then 1678 msg='wrong sizes (pawcprj_reorder)!' 1679 MSG_BUG(msg) 1680 end if 1681 1682 !Nothing to do when the atoms are already sorted 1683 iexit=1;ii=0 1684 do while (iexit==1.and.ii<n1atindx) 1685 ii=ii+1 1686 if (atm_indx(ii)/=ii) iexit=0 1687 end do 1688 if (iexit==1) return 1689 1690 LIBPAW_ALLOCATE(nlmn,(n1cprj)) 1691 do ii=1,n1cprj 1692 nlmn(ii)=cprj(ii,1)%nlmn 1693 end do 1694 ncpgr=cprj(1,1)%ncpgr 1695 LIBPAW_DATATYPE_ALLOCATE(cprj_tmp,(n1cprj,n2cprj)) 1696 call pawcprj_alloc(cprj_tmp,ncpgr,nlmn) 1697 call pawcprj_copy(cprj,cprj_tmp) 1698 call pawcprj_free(cprj) 1699 1700 do jj=1,n2cprj 1701 do ii=1,n1cprj 1702 kk=atm_indx(ii) 1703 cprj(kk,jj)%nlmn=nlmn(ii) 1704 cprj(kk,jj)%ncpgr=ncpgr 1705 LIBPAW_ALLOCATE(cprj(kk,jj)%cp,(2,nlmn(ii))) 1706 cprj(kk,jj)%cp(:,:)=cprj_tmp(ii,jj)%cp(:,:) 1707 if (ncpgr>0) then 1708 LIBPAW_ALLOCATE(cprj(kk,jj)%dcp,(2,ncpgr,nlmn(ii))) 1709 cprj(kk,jj)%dcp(:,:,:)=cprj_tmp(ii,jj)%dcp(:,:,:) 1710 end if 1711 end do 1712 end do 1713 1714 call pawcprj_free(cprj_tmp) 1715 LIBPAW_DATATYPE_DEALLOCATE(cprj_tmp) 1716 LIBPAW_DEALLOCATE(nlmn) 1717 1718 end subroutine pawcprj_reorder
m_pawcprj/pawcprj_set_zero [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_set_zero
FUNCTION
Set to zero all arrays in a cprj datastructure
SIDE EFFECTS
cprj(:,:) <type(pawcprj_type)>= cprj datastructure
PARENTS
ctocprj,dfpt_cgwf,m_fock
CHILDREN
xmpi_sum
SOURCE
286 subroutine pawcprj_set_zero(cprj) 287 288 289 !This section has been created automatically by the script Abilint (TD). 290 !Do not modify the following lines by hand. 291 #undef ABI_FUNC 292 #define ABI_FUNC 'pawcprj_set_zero' 293 !End of the abilint section 294 295 implicit none 296 297 !Arguments ------------------------------------ 298 !scalars 299 !arrays 300 type(pawcprj_type),intent(inout) :: cprj(:,:) 301 302 !Local variables------------------------------- 303 !scalars 304 integer :: ii,jj,n1dim,n2dim 305 306 ! ************************************************************************* 307 308 n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2) 309 310 do jj=1,n2dim 311 do ii=1,n1dim 312 if (cprj(ii,jj)%nlmn>0) cprj(ii,jj)%cp(:,:)=zero 313 if (cprj(ii,jj)%ncpgr>0) cprj(ii,jj)%dcp(:,:,:)=zero 314 end do 315 end do 316 317 end subroutine pawcprj_set_zero
m_pawcprj/pawcprj_symkn [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_symkn
FUNCTION
compute cprj for a given band and k point based on cprj at a symmetry-related k point.
INPUTS
cprj_ikn (pawcprj_type) :: cprj for a single band and k point, typically a k point in the IBZ cprj_sym(4,nsym,natom) :: 1:3 shift, and 4 final atom, of symmetry isym operating on iatom (S^{-1}(R - t) = r0 + L, see symatm.F90 dimlmn(natom) :: ln dimension of each atom iband :: number of bands to treat, use -1 to treat all nband bands indlmn(6,lmnmax,ntypat) :: n,l,m dimensions for each atom type (see psps type) isym :: symmetry element used in current application itim :: 1 if time reversal also used, 0 else kpt(3) :: kpt vector used lmax :: max l value lmnmax :: max lmn value mband :: maximum number of bands natom :: number of atoms in cell nband :: number of bands in cprj_ikn nspinor :: number of spinors nsym :: total number of symmetry elements ntypat :: number of types of atoms typat(natom) :: type of each atom zarot(2*lmax+1,2*lmax+1,lmax+1,nsym) :: elements of rotation matrix for angular momentum states and symmetry operations. See m_paw_sphharm/setsym_ylm.
OUTPUT
cprj_fkn (pawcprj_type) :: cprj for a single band and k point where the k point is related to the input k point by a symmetry operation
SIDE EFFECTS
NOTES
This routine is based on M. Giantomassi's doctoral dissertation, formula 7.77. It is not clear whether it is implemented correctly for nonsymmorphic symmetries.
PARENTS
berryphase_new,cgwf,m_fock,make_grad_berry
CHILDREN
xmpi_sum
SOURCE
779 subroutine pawcprj_symkn(cprj_fkn,cprj_ikn,cprj_sym,dimlmn,iband,indlmn,& 780 & isym,itim,kpt,lmax,lmnmax,mband,natom,nband,nspinor,nsym,ntypat,& 781 & typat,zarot) 782 783 784 !This section has been created automatically by the script Abilint (TD). 785 !Do not modify the following lines by hand. 786 #undef ABI_FUNC 787 #define ABI_FUNC 'pawcprj_symkn' 788 !End of the abilint section 789 790 implicit none 791 792 !Arguments--------------------------- 793 !scalars 794 integer,intent(in) :: iband,isym,itim,lmax,lmnmax,mband 795 integer,intent(in) :: natom,nband,nspinor,nsym,ntypat 796 797 !arrays 798 integer,intent(in) :: cprj_sym(4,nsym,natom),dimlmn(natom) 799 integer,intent(in) :: indlmn(6,lmnmax,ntypat),typat(natom) 800 real(dp),intent(in) :: kpt(3) 801 real(dp),intent(in) :: zarot(2*lmax+1,2*lmax+1,lmax+1,nsym) 802 type(pawcprj_type),intent(in) :: cprj_ikn(natom,mband*nspinor) 803 type(pawcprj_type),intent(inout) :: cprj_fkn(natom,mband*nspinor) !vz_i 804 805 !Local variables--------------------------- 806 !scalars 807 integer :: iatm,iatom, ibct, ibnd, ibsp, ibst, icpgr, iin, il, il0, im 808 integer :: ilmn, iln, iln0, ilpm, indexi, ispinor, itypat, jatm,jatom, mm, nlmn 809 real(dp) :: kdotL, phr, phi 810 !arrays 811 real(dp) :: rl(3), t1(2), t2(2) 812 813 ! ************************************************************************* 814 815 if (iband == -1) then 816 ibst = 1 817 ibnd = nband 818 else 819 ibst = iband 820 ibnd = iband 821 end if 822 823 do iatom = 1, natom 824 iatm=iatom 825 itypat = typat(iatom) 826 nlmn = dimlmn(iatm) 827 jatom = cprj_sym(4,isym,iatom) 828 jatm=jatom 829 rl(:) = cprj_sym(1:3,isym,iatom) 830 kdotL = dot_product(rl,kpt) 831 phr = cos(two_pi*kdotL) 832 phi = sin(two_pi*kdotL) 833 834 il0 = -1; iln0 = -1; indexi = 1 835 do ilmn = 1, nlmn 836 837 il = indlmn(1,ilmn,itypat) 838 im = indlmn(2,ilmn,itypat) 839 iin = indlmn(3,ilmn,itypat) 840 iln = indlmn(5,ilmn,itypat) 841 ilpm = 1 + il + im 842 if (iln /= iln0) indexi = indexi + 2*il0 + 1 843 844 do ibct = ibst, ibnd 845 846 do ispinor = 1, nspinor 847 848 ibsp = nspinor*(ibct-1) + ispinor 849 850 t1(:) = zero 851 do mm = 1, 2*il+1 852 t1(1) = t1(1) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%cp(1,indexi+mm) 853 t1(2) = t1(2) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%cp(2,indexi+mm) 854 end do 855 t2(1) = t1(1)*phr - t1(2)*phi 856 t2(2) = t1(2)*phr + t1(1)*phi 857 858 if (itim == 1) t2(2) = -t2(2) 859 860 cprj_fkn(iatm,ibsp)%cp(1,ilmn) = t2(1) 861 cprj_fkn(iatm,ibsp)%cp(2,ilmn) = t2(2) 862 863 ! do same transformations for gradients of cprj_ikn 864 ! note that ncpgr = 0 if no gradients present so this loop will not be executed 865 ! in this case 866 867 do icpgr = 1, cprj_ikn(jatom,ibsp)%ncpgr 868 t1(:) = zero 869 870 do mm = 1, 2*il+1 871 t1(1) = t1(1) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%dcp(1,icpgr,indexi+mm) 872 t1(2) = t1(2) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%dcp(2,icpgr,indexi+mm) 873 end do 874 875 t2(1) = t1(1)*phr - t1(2)*phi 876 t2(2) = t1(2)*phr + t1(1)*phi 877 878 if (itim == 1) t2(2) = -t2(2) 879 880 cprj_fkn(iatm,ibsp)%dcp(1,icpgr,ilmn) = t2(1) 881 cprj_fkn(iatm,ibsp)%dcp(2,icpgr,ilmn) = t2(2) 882 883 end do ! end loop over ncpgr 884 885 end do ! end loop over nspinor 886 887 end do ! end loop over bands 888 889 il0 = il; iln0 = iln 890 end do ! end loop over ilmn 891 end do ! end loop over atoms 892 893 end subroutine pawcprj_symkn
m_pawcprj/pawcprj_transpose [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_transpose
FUNCTION
Transpose a cprj datastructure FOR A GIVEN (K,SPIN) in order to change the parallel distribution from atom to band (or the contrary). At input, cprj is distributed over bands (or atoms); at output, it is distributed over atoms (or bands)
INPUTS
cprjin(n1indim,n2indim)<pawcprj_type>=the input cprj datastructure cprj_bandpp=number of bands to be treated simultaneoulsy by a processor natom=number of atoms in cell nband=number of bands nspinor=number of spinorial components spaceComm=MPI Communicator.
OUTPUT
cprjout(n1outdim,n2outdim)<pawcprj_type>=the output cprj datastructure with another distribution
NOTES
On the dimensions: To transfer cprj from band distribution to atom distribution, dimensions should be: n1indim =natom n2indim =nband/nproc*nspinor n1outdim=natom/nproc n2outdim=nband*nspinor To transfer cprj from atom distribution to band distribution, dimensions should be: n1indim =natom n2indim =nband/nproc*nspinor n1outdim=natom/nproc n2outdim=nband*nspinor
PARENTS
CHILDREN
xmpi_sum
SOURCE
2503 subroutine pawcprj_transpose(cprjin,cprjout,cprj_bandpp,natom,nband,nspinor,spaceComm) 2504 2505 2506 !This section has been created automatically by the script Abilint (TD). 2507 !Do not modify the following lines by hand. 2508 #undef ABI_FUNC 2509 #define ABI_FUNC 'pawcprj_transpose' 2510 !End of the abilint section 2511 2512 implicit none 2513 2514 !Arguments------------------------------------- 2515 !scalars 2516 integer :: cprj_bandpp,natom,nband,nspinor,spaceComm 2517 !arrays 2518 type(pawcprj_type),intent(in) :: cprjin(:,:) 2519 type(pawcprj_type),intent(out) :: cprjout(:,:) 2520 2521 !Local variables------------------------------- 2522 !scalars 2523 integer :: bpp,buf_indx 2524 integer :: iashft,iatom,iatom_max_sd,iatom_max_rc,iatom_1,iatom_2,iatm1_sd,iatm1_rc,iatm2_sd,iatm2_rc 2525 integer :: ib,iband,iband_1,iband_2,iband_shift,iblock_atom,iblock_band,ibshft 2526 integer :: ierr,ip,ispinor,me,nba,nbb,nbnp_sd,nbnp_rc,ncpgr,nlmn,np 2527 integer :: rbufsize,sbufsize,size11,size12,size21,size22,transpose_mode 2528 character(len=100) :: msg 2529 !arrays 2530 integer,allocatable :: cprjsz_atom(:),cprjsz_block(:,:) 2531 integer,allocatable,target :: count_atom(:),count_band(:),displ_atom(:),displ_band(:) 2532 integer,pointer :: scount(:),sdispl(:),rcount(:),rdispl(:) 2533 real(dp),allocatable :: rbuf(:),sbuf(:) 2534 2535 ! ************************************************************************* 2536 2537 !MPI data 2538 me = xmpi_comm_rank(spaceComm) 2539 np = xmpi_comm_size(spaceComm) 2540 2541 !Nothing to do if nprocs=1 2542 if (np==1) then 2543 call pawcprj_copy(cprjin,cprjout) 2544 return 2545 end if 2546 2547 !Compute bloc sizes 2548 bpp=cprj_bandpp 2549 nba=natom/np;if (mod(natom,np)/=0) nba=nba+1 2550 nbb=nband/(np*bpp) 2551 2552 !Check sizes, select direction of transposition 2553 transpose_mode=0 2554 size11=size(cprjin,1);size12=size(cprjin,2) 2555 size21=size(cprjout,1);size22=size(cprjout,2) 2556 if (size11==natom.and.size12==nbb*bpp*nspinor.and.& 2557 & size21==nba.and.size22==nband*nspinor) then 2558 transpose_mode=1 2559 else if (size11==nba.and.size12==nband*nspinor.and.& 2560 & size21==natom.and.size22==nbb*bpp*nspinor) then 2561 else 2562 msg='wrong cprjin/cprjout sizes (pawcprj_transpose)!' 2563 MSG_BUG(msg) 2564 end if 2565 2566 !Compute size of atom bloc (wr to cprj) 2567 LIBPAW_ALLOCATE(cprjsz_atom,(natom)) 2568 LIBPAW_ALLOCATE(cprjsz_block,(np,nba)) 2569 cprjsz_atom=0;cprjsz_block=0 2570 if (transpose_mode==1) then 2571 do iatom=1,natom 2572 cprjsz_atom(iatom)=2*cprjin(iatom,1)%nlmn*(1+cprjin(iatom,1)%ncpgr) 2573 end do 2574 else 2575 do iblock_atom=1,nba 2576 iatom=(iblock_atom-1)*np+1+me 2577 if (iatom<=natom) cprjsz_atom(iatom)=2*cprjin(iblock_atom,1)%nlmn*(1+cprjin(iblock_atom,1)%ncpgr) 2578 end do 2579 call xmpi_sum(cprjsz_atom,spaceComm,ierr) 2580 end if 2581 do iblock_atom=1,nba 2582 iashft=(iblock_atom-1)*np 2583 iatom_1=iashft+1;iatom_2=iashft+np 2584 if (iatom_1>natom) cycle 2585 if (iatom_2>natom) iatom_2=natom 2586 do iatom=iatom_1,iatom_2 2587 cprjsz_block(iatom-iashft,iblock_atom)=cprjsz_atom(iatom)+2 ! +2 for nlmn et ncpgr 2588 end do 2589 end do 2590 LIBPAW_DEALLOCATE(cprjsz_atom) 2591 2592 !Allocations for MPI_ALLTOALL 2593 LIBPAW_ALLOCATE(count_atom,(np)) 2594 LIBPAW_ALLOCATE(displ_atom,(np)) 2595 LIBPAW_ALLOCATE(count_band,(np)) 2596 LIBPAW_ALLOCATE(displ_band,(np)) 2597 2598 !Loop on blocks of bands 2599 do iblock_band=1,nbb !(note: np divides nband) 2600 ibshft=(iblock_band-1)*np*bpp 2601 iband_1=ibshft+1;iband_2=ibshft+np*bpp 2602 if (iband_1>nband.or.iband_2>nband) cycle ! for security 2603 2604 ! Loop on blocks of atoms 2605 do iblock_atom=1,nba 2606 iashft=(iblock_atom-1)*np 2607 iatom_1=iashft+1;iatom_2=iashft+np 2608 if (iatom_1>natom) cycle 2609 if (iatom_2>natom) iatom_2=natom 2610 2611 ! Computation of displacements and sizes of blocks when data are band-distributed 2612 count_band(1)=cprjsz_block(1,iblock_atom)*nspinor*bpp;displ_band(1)=0 2613 do ip=2,np 2614 count_band(ip)=cprjsz_block(ip,iblock_atom)*nspinor*bpp 2615 displ_band(ip)=displ_band(ip-1)+count_band(ip-1) 2616 end do 2617 2618 ! Computation of displacements and sizes of blocks when data are atom-distributed 2619 count_atom(1)=cprjsz_block(1+me,iblock_atom)*bpp*nspinor;displ_atom(1)=0 2620 do ip=2,np 2621 count_atom(ip)=count_atom(1) 2622 displ_atom(ip)=displ_atom(ip-1)+count_atom(ip-1) 2623 end do 2624 2625 ! According to transposition mode, select 2626 ! - displacements and sizes of blocks 2627 ! - shifts in arrays 2628 if (transpose_mode==1) then 2629 scount => count_band ; sdispl => displ_band 2630 rcount => count_atom ; rdispl => displ_atom 2631 nbnp_sd=bpp;nbnp_rc=np*bpp 2632 iatm1_sd=iatom_1;iatm2_sd=iatom_2 2633 iatm1_rc=iblock_atom;iatm2_rc=iatm1_rc 2634 iatom_max_sd=iatom_2;iatom_max_rc=iashft+1+me 2635 else 2636 scount => count_atom ; sdispl => displ_atom 2637 rcount => count_band ; rdispl => displ_band 2638 nbnp_sd=np*bpp;nbnp_rc=bpp 2639 iatm1_sd=iblock_atom;iatm2_sd=iatm1_sd 2640 iatm1_rc=iatom_1;iatm2_rc=iatom_2 2641 iatom_max_sd=iashft+1+me;iatom_max_rc=iatom_2 2642 end if 2643 2644 ! Allocation of buffers 2645 sbufsize=sdispl(np)+scount(np) 2646 rbufsize=rdispl(np)+rcount(np) 2647 LIBPAW_ALLOCATE(sbuf,(sbufsize)) 2648 LIBPAW_ALLOCATE(rbuf,(rbufsize)) 2649 2650 ! Coying of input cprj to buffer for sending 2651 buf_indx=0 2652 iband_shift=(iblock_band-1)*nbnp_sd-1 2653 if (iatom_max_sd<=natom) then 2654 do iatom=iatm1_sd,iatm2_sd 2655 do ib=1,nbnp_sd 2656 iband=(iband_shift+ib)*nspinor 2657 do ispinor=1,nspinor 2658 iband=iband+1 2659 nlmn=cprjin(iatom,iband)%nlmn;ncpgr=cprjin(iatom,iband)%ncpgr 2660 sbuf(buf_indx+1)=dble(nlmn) ;buf_indx=buf_indx+1 2661 sbuf(buf_indx+1)=dble(ncpgr);buf_indx=buf_indx+1 2662 sbuf(buf_indx+1:buf_indx+2*nlmn)=reshape(cprjin(iatom,iband)%cp(1:2,1:nlmn),(/2*nlmn/)) 2663 buf_indx=buf_indx+2*nlmn 2664 if (ncpgr>0) then 2665 sbuf(buf_indx+1:buf_indx+2*ncpgr*nlmn)=reshape(cprjin(iatom,iband)%dcp(1:2,1:ncpgr,1:nlmn),(/2*ncpgr*nlmn/)) 2666 buf_indx=buf_indx+2*ncpgr*nlmn 2667 end if 2668 end do 2669 end do 2670 end do 2671 end if 2672 if (buf_indx/=sbufsize) then 2673 msg='wrong buffer size for sending (pawcprj_transpose)!' 2674 MSG_BUG(msg) 2675 end if 2676 2677 ! Main call to MPI_ALLTOALL 2678 call xmpi_alltoallv(sbuf,scount,sdispl,rbuf,rcount,rdispl,spaceComm,ierr) 2679 2680 ! Retrieving of output cprj for received buffer 2681 buf_indx=0 2682 iband_shift=(iblock_band-1)*nbnp_rc-1 2683 if (iatom_max_rc<=natom) then 2684 do iatom=iatm1_rc,iatm2_rc 2685 do ib=1,nbnp_rc 2686 iband=(iband_shift+ib)*nspinor 2687 do ispinor=1,nspinor 2688 iband=iband+1 2689 nlmn =int(rbuf(buf_indx+1));buf_indx=buf_indx+1 2690 ncpgr=int(rbuf(buf_indx+1));buf_indx=buf_indx+1 2691 cprjout(iatom,iband)%nlmn=nlmn;cprjout(iatom,iband)%ncpgr=ncpgr 2692 cprjout(iatom,iband)%cp(1:2,1:nlmn)=reshape(rbuf(buf_indx+1:buf_indx+2*nlmn),(/2,nlmn/)) 2693 buf_indx=buf_indx+2*nlmn 2694 if (ncpgr>0) then 2695 cprjout(iatom,iband)%dcp(1:2,1:ncpgr,1:nlmn)=reshape(rbuf(buf_indx+1:buf_indx+2*nlmn*ncpgr),(/2,ncpgr,nlmn/)) 2696 buf_indx=buf_indx+2*nlmn*ncpgr 2697 end if 2698 end do 2699 end do 2700 end do 2701 else 2702 cprjout(iatom,iband)%nlmn=0;cprjout(iatom,iband)%ncpgr=0 2703 end if 2704 if (buf_indx/=rbufsize) then 2705 msg='wrong buffer size for receiving (pawcprj_transpose)!' 2706 MSG_BUG(msg) 2707 end if 2708 2709 ! Deallocation of buffers 2710 LIBPAW_DEALLOCATE(sbuf) 2711 LIBPAW_DEALLOCATE(rbuf) 2712 2713 ! End of loops 2714 end do ! do iblock_atom 2715 end do ! do iblock_atom 2716 2717 !Free memory 2718 LIBPAW_DEALLOCATE(count_atom) 2719 LIBPAW_DEALLOCATE(displ_atom) 2720 LIBPAW_DEALLOCATE(count_band) 2721 LIBPAW_DEALLOCATE(displ_band) 2722 LIBPAW_DEALLOCATE(cprjsz_block) 2723 nullify(scount,rcount,sdispl,rdispl) 2724 2725 end subroutine pawcprj_transpose
m_pawcprj/pawcprj_type [ Types ]
[ Top ] [ m_pawcprj ] [ Types ]
NAME
pawcprj_type
FUNCTION
This structured datatype contains <p_lmn|Cnk> projected scalars and derivatives where |p_lmn> are non-local projectors for a given atom |Cnk> is a wave function Used only for PAW calculations.
SOURCE
53 type,public :: pawcprj_type 54 55 !Integer scalars 56 57 integer :: ncpgr=0 58 ! Number of gradients of cp=<p_lmn|Cnk> 59 60 integer :: nlmn=0 61 ! Number of (l,m,n) non-local projectors 62 63 !Real (real(dp)) arrays 64 65 real(dp), allocatable :: cp (:,:) 66 ! cp(2,nlmn) 67 ! <p_lmn|Cnk> projected scalars for a given atom and wave function 68 69 real(dp), allocatable :: dcp (:,:,:) 70 ! dcp(2,ncpgr,nlmn) 71 ! derivatives of <p_lmn|Cnk> projected scalars for a given atom and wave function 72 73 end type pawcprj_type 74 75 !public procedures. 76 public :: pawcprj_alloc ! Allocation 77 public :: pawcprj_free ! Deallocation 78 public :: pawcprj_set_zero ! Set to zero all arrays in a cprj datastructure 79 public :: pawcprj_copy ! Copy a cprj datastructure into another 80 public :: pawcprj_axpby ! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:) 81 public :: pawcprj_zaxpby ! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:), alpha and beta are COMPLEX scalars 82 public :: pawcprj_conjg ! cprj(:,:) <- conjugate(cprj(:,:)) 83 public :: pawcprj_symkn ! construct cprj from that at a symmetry related k point 84 public :: pawcprj_lincom ! Compute a LINear COMbination of cprj datastructure: 85 public :: pawcprj_output ! Output a cprj. Useful for debugging. 86 public :: pawcprj_get ! Read the cprj for a given k-point from memory or from a temporary file 87 public :: pawcprj_put ! Write the cprj for a given set of (n,k) into memory or into a temporary file 88 public :: pawcprj_reorder ! Change the order of a cprj datastructure 89 public :: pawcprj_mpi_allgather ! Perform MPI_ALLGATHER on a pawcprj_type inside a MPI communicator. 90 public :: pawcprj_bcast ! Broadcast a pawcprj_type from master to all nodes inside a MPI communicator. 91 public :: pawcprj_transpose ! Transpose a cprj datastructure FOR A GIVEN (K,SPIN) 92 public :: pawcprj_gather_spin ! Collect spin distributed cprjs. 93 public :: pawcprj_mpi_exch ! Exchange a pawcprj_type between two processors inside a MPI communicator. 94 public :: pawcprj_mpi_send ! Send a pawcprj_type inside a MPI communicator. 95 public :: pawcprj_mpi_recv ! Receive a pawcprj_type inside a MPI communicator. 96 public :: pawcprj_mpi_sum ! Perform MPI_SUM on a pawcprj_type inside a MPI communicator. 97 public :: pawcprj_getdim ! Returns the number of lmn components in the <p_{lmn}^i|\psi> for the i-th atom. 98 public :: paw_overlap ! Compute the onsite contribution to the overlap between two states.
m_pawcprj/pawcprj_zaxpby [ Functions ]
[ Top ] [ m_pawcprj ] [ Functions ]
NAME
pawcprj_zaxpby
FUNCTION
Apply ZAXPBY (blas-like) operation with 2 cprj datastructures: cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:) alpha and beta are COMPLEX scalars
INPUTS
alpha(2),beta(2)= alpha,beta COMPLEX factors cprjx(:,:) <type(pawcprj_type)>= input cprjx datastructure
SIDE EFFECTS
cprjy(:,:) <type(pawcprj_type)>= input/output cprjy datastructure
PARENTS
corrmetalwf1,extrapwf
CHILDREN
xmpi_sum
SOURCE
601 subroutine pawcprj_zaxpby(alpha,beta,cprjx,cprjy) 602 603 604 !This section has been created automatically by the script Abilint (TD). 605 !Do not modify the following lines by hand. 606 #undef ABI_FUNC 607 #define ABI_FUNC 'pawcprj_zaxpby' 608 !End of the abilint section 609 610 implicit none 611 612 !Arguments ------------------------------------ 613 !scalars 614 real(dp),intent(in) :: alpha(2),beta(2) 615 !arrays 616 type(pawcprj_type),intent(in) :: cprjx(:,:) 617 type(pawcprj_type),intent(inout) :: cprjy(:,:) 618 619 !Local variables------------------------------- 620 !scalars 621 integer :: ii,jj,kk,ll,n1dimx,n1dimy,n2dimx,n2dimy,ncpgrx,ncpgry,nlmn 622 real(dp) :: cp1,cp2,norma,normb 623 character(len=500) :: msg 624 625 ! ************************************************************************* 626 627 norma=alpha(1)**2+alpha(2)**2 628 normb=beta(1) **2+beta(2) **2 629 n1dimy=size(cprjy,dim=1);n2dimy=size(cprjy,dim=2);ncpgry=cprjy(1,1)%ncpgr 630 if (norma>tol16) then 631 n1dimx=size(cprjx,dim=1);n2dimx=size(cprjx,dim=2);ncpgrx=cprjx(1,1)%ncpgr 632 msg = "" 633 if (n1dimx/=n1dimy) msg = TRIM(msg)//"Error in pawcprj_zaxpby: n1 wrong sizes !"//ch10 634 if (n2dimx/=n2dimy) msg = TRIM(msg)//"Error in pawcprj_zaxpby: n2 wrong sizes !"//ch10 635 if (ncpgrx/=ncpgry) msg = TRIM(msg)//"Error in pawcprj_zaxpby: ncpgr wrong sizes !"//ch10 636 if (LEN_TRIM(msg) > 0) then 637 MSG_ERROR(msg) 638 end if 639 end if 640 641 if (norma<=tol16) then 642 do jj=1,n2dimy 643 do ii=1,n1dimy 644 nlmn=cprjy(ii,jj)%nlmn 645 do kk=1,nlmn 646 cp1=beta(1)*cprjy(ii,jj)%cp(1,kk)-beta(2)*cprjy(ii,jj)%cp(2,kk) 647 cp2=beta(1)*cprjy(ii,jj)%cp(2,kk)+beta(2)*cprjy(ii,jj)%cp(1,kk) 648 cprjy(ii,jj)%cp(1,kk)=cp1 649 cprjy(ii,jj)%cp(2,kk)=cp2 650 end do 651 end do 652 end do 653 if (ncpgry>0) then 654 do jj=1,n2dimy 655 do ii=1,n1dimy 656 nlmn=cprjy(ii,jj)%nlmn 657 do kk=1,nlmn 658 do ll=1,ncpgry 659 cp1=beta(1)*cprjy(ii,jj)%dcp(1,ll,kk)-beta(2)*cprjy(ii,jj)%dcp(2,ll,kk) 660 cp2=beta(1)*cprjy(ii,jj)%dcp(2,ll,kk)+beta(2)*cprjy(ii,jj)%dcp(1,ll,kk) 661 cprjy(ii,jj)%dcp(1,ll,kk)=cp1 662 cprjy(ii,jj)%dcp(2,ll,kk)=cp2 663 end do 664 end do 665 end do 666 end do 667 end if 668 else if (normb<=tol16) then 669 do jj=1,n2dimx 670 do ii=1,n1dimx 671 nlmn=cprjx(ii,jj)%nlmn 672 cprjy(ii,jj)%nlmn=nlmn 673 do kk=1,nlmn 674 cprjy(ii,jj)%cp(1,kk)=alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk) 675 cprjy(ii,jj)%cp(2,kk)=alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk) 676 end do 677 end do 678 end do 679 if (ncpgrx>0) then 680 do jj=1,n2dimx 681 do ii=1,n1dimx 682 nlmn=cprjx(ii,jj)%nlmn 683 do kk=1,nlmn 684 cprjy(ii,jj)%dcp(1,1:ncpgrx,kk)=alpha(1)*cprjx(ii,jj)%dcp(1,1:ncpgrx,kk) & 685 & -alpha(2)*cprjx(ii,jj)%dcp(2,1:ncpgrx,kk) 686 cprjy(ii,jj)%dcp(2,1:ncpgrx,kk)=alpha(1)*cprjx(ii,jj)%dcp(2,1:ncpgrx,kk) & 687 & +alpha(2)*cprjx(ii,jj)%dcp(1,1:ncpgrx,kk) 688 end do 689 end do 690 end do 691 end if 692 else 693 do jj=1,n2dimx 694 do ii=1,n1dimx 695 nlmn=cprjx(ii,jj)%nlmn 696 cprjy(ii,jj)%nlmn =nlmn 697 do kk=1,nlmn 698 cp1=alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk) & 699 & +beta(1) *cprjy(ii,jj)%cp(1,kk)-beta(2) *cprjy(ii,jj)%cp(2,kk) 700 cp2=alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk) & 701 & +beta(1) *cprjy(ii,jj)%cp(2,kk)+beta(2) *cprjy(ii,jj)%cp(1,kk) 702 cprjy(ii,jj)%cp(1,kk)=cp1 703 cprjy(ii,jj)%cp(2,kk)=cp2 704 end do 705 end do 706 end do 707 if (ncpgrx>0) then 708 do jj=1,n2dimx 709 do ii=1,n1dimx 710 nlmn=cprjx(ii,jj)%nlmn 711 do kk=1,nlmn 712 do ll=1,ncpgrx 713 cp1=alpha(1)*cprjx(ii,jj)%dcp(1,ll,kk)-alpha(2)*cprjx(ii,jj)%dcp(2,ll,kk) & 714 & +beta(1) *cprjy(ii,jj)%dcp(1,ll,kk)-beta(2) *cprjy(ii,jj)%dcp(2,ll,kk) 715 cp2=alpha(1)*cprjx(ii,jj)%dcp(2,ll,kk)+alpha(2)*cprjx(ii,jj)%dcp(1,ll,kk) & 716 & +beta(1) *cprjy(ii,jj)%dcp(2,ll,kk)+beta(2) *cprjy(ii,jj)%dcp(1,ll,kk) 717 cprjy(ii,jj)%dcp(1,ll,kk)=cp1 718 cprjy(ii,jj)%dcp(2,ll,kk)=cp2 719 end do 720 end do 721 end do 722 end do 723 end if 724 end if 725 726 end subroutine pawcprj_zaxpby