TABLE OF CONTENTS
ABINIT/cprj_rotate [ Functions ]
NAME
cprj_rotate
FUNCTION
Compute cprj_nk = \sum_m z_m cprj_mk where z_m is an array of complex values. The input is overwritten.
INPUTS
SIDE EFFECTS
SOURCE
1003 subroutine cprj_rotate(cprj_in,evec,dimcprj,natom,nband,nspinor) 1004 1005 !Arguments ------------------------------- 1006 !scalars 1007 integer,intent(in) :: natom,nband,nspinor 1008 !arrays 1009 integer,intent(in) :: dimcprj(:) 1010 real(dp) :: evec(:,:) 1011 type(pawcprj_type),intent(inout) :: cprj_in(natom,nspinor*nband) 1012 1013 !Local variables------------------------------- 1014 !scalars 1015 integer :: iband,ncpgr 1016 !arrays 1017 ! real(dp) :: tsec(2) 1018 real(dp) :: z_tmp(2,nband) 1019 type(pawcprj_type),pointer :: cprj_iband(:,:) 1020 type(pawcprj_type),allocatable,target :: cprj_tmp(:,:) 1021 1022 1023 ! ********************************************************************* 1024 1025 DBG_ENTER('COLL') 1026 1027 ! call timab(1211,1,tsec) 1028 1029 ncpgr=cprj_in(1,1)%ncpgr 1030 ABI_MALLOC(cprj_tmp,(natom,nspinor*nband)) 1031 call pawcprj_alloc(cprj_tmp,ncpgr,dimcprj) 1032 1033 do iband=1,nband 1034 z_tmp = reshape(evec(:,iband),(/2,nband/)) 1035 cprj_iband => cprj_tmp(:,nspinor*(iband-1)+1:nspinor*iband) 1036 call pawcprj_lincom(z_tmp,cprj_in,cprj_iband,nband) 1037 end do 1038 1039 call pawcprj_copy(cprj_tmp,cprj_in) 1040 call pawcprj_free(cprj_tmp) 1041 ABI_FREE(cprj_tmp) 1042 1043 ! call timab(1211,2,tsec) 1044 1045 DBG_EXIT('COLL') 1046 1047 end subroutine cprj_rotate
ABINIT/ctocprj [ Functions ]
NAME
ctocprj
FUNCTION
Compute all <Proj_i|Cnk> for every wave function |Cnk> expressed in reciprocal space. |Proj_i> are non-local projectors (for each atom and each l,m,n) Can also compute derivatives of <Proj_i|Cnk> wrt to several parameters
INPUTS
atindx(natom)=index table for atoms cg(2,mcg)=planewave coefficients of wavefunctions choice: chooses derivatives to compute: =1 => no derivatives =2 => 1st derivatives with respect to atomic position(s) =3 => 1st derivatives with respect to strain(s) =23=> 1st derivatives with respect to strain(s) and atm pos =4 => 2nd derivatives with respect to atomic pos. =24=> 1st and 2nd derivatives with respect to atomic pos. =5 => derivatives with respect to k wavevector =6 => 2nd derivatives with respect to strain and atm. pos. gmet(3,3)=reciprocal space metric tensor in bohr**-2 gprimd(3,3)=dimensional reciprocal space primitive translations iatom= if <=0, cprj=<p_i|Cnk> are computed for all atoms 1...natom if >0 cprj=<p_i|Cnk> are computed only for atom with index iatom idir=direction of the derivative, i.e. dir. of - atom to be moved in the case choice=2 - strain component in the case choice=3 - k point direction in the case choice=5 Compatible only with choice=2,3,5; if idir=0, all derivatives are computed iorder_cprj=0 if output cprj=<p_i|Cnk> are sorted by atom type (first all elements of atom type 1, followed by those of atom type 2 and so on). 1 if output cprj=<p_i|Cnk> are sorted according to the variable typat in the main input file istwfk(nkpt)=option parameter that describes the storage of wfs kg(3,mpw*mkmem)=reduced planewave coordinates kpt(3,nkpt)=reduced coordinates of k points. mcg=size of wave-functions array (cg) =mpw*nspinor*mband_mem*mkmem*nsppol mcprj=size of projected wave-functions array (cprj) =nspinor*mband*mkmem*nsppol mgfft=maximum size of 1D FFTs mkmem=number of k points treated by this node. mpi_enreg=information about MPI parallelization mpsang=1+maximum angular momentum for nonlocal pseudopotentials mpw=maximum dimensioned size of npw my_nsppol=number of spin components in memory for current MPI process natom=number of atoms in cell nattyp(ntypat)= # atoms of each type nband(nkpt*nsppol)=number of bands at this k point for that spin polarization mband_mem=max number of bands for this processor (in case of band parallelism) ncprj=1st dim. of cprj array (natom if iatom<=0, 1 if iatom>0) ngfft(18)=contain all needed information about 3D FFT, see ~ABINIT/Infos/vargs.htm#ngfft nkpt=number of k points nloalg(3)=governs the choice of the algorithm for nonlocal operator npwarr(nkpt)=number of planewaves in basis at this k point nspinor=number of spinorial components of the wavefunctions nsppol=1 for unpolarized, 2 for spin-polarized ntypat=number of types of atoms in unit cell paral_kgb= 1 if kpt-band-FFT is activated ph1d(2,3*(2*mgfft+1)*natom)=1-dim structure factor phase information psps <type(pseudopotential_type)>=variables related to pseudopotentials rmet(3,3)=real space metric (bohr**2) tim_ctocprj=timing code of the calling routine typat(natom)= types of atoms uncp=unit number for <P_lmn|Cnk> data (if used) xred(3,natom)=reduced dimensionless atomic coordinates ylm(mpw*mkmem,mpsang*mpsang)=real spherical harmonics for each G and k point ylmgr(npw*mkmem,nylmgr,mpsang*mpsang*useylmgr)=gradients of real spherical harmonics wrt (k+G)
OUTPUT
cprj(ncprj,mcprj) <type(pawcprj_type)>= projected input wave functions <Proj_i|Cnk> with NL projectors Usually ncprj=natom
SOURCE
447 subroutine ctocprj(atindx,cg,choice,cprj,gmet,gprimd,iatom,idir,& 448 & iorder_cprj,istwfk,kg,kpt,mcg,mcprj,mgfft,mkmem,mpi_enreg,mpsang,& 449 & mpw,natom,nattyp,nband,ncprj,ngfft,nkpt,nloalg,npwarr,nspinor,& 450 & nsppol,my_nsppol,ntypat,paral_kgb,ph1d,psps,rmet,typat,ucvol,uncp,xred,ylm,ylmgr) 451 452 !Arguments ------------------------------- 453 !scalars 454 integer,intent(in) :: choice,iatom,idir,iorder_cprj,mcg,mcprj,mgfft,mkmem,mpsang,mpw 455 integer,intent(in) :: my_nsppol,natom,ncprj,nkpt,nspinor,nsppol,ntypat,paral_kgb,uncp 456 !TODO : distribute cprj over bands as well 457 real(dp),intent(in) :: ucvol 458 type(MPI_type),intent(in) :: mpi_enreg 459 type(pseudopotential_type),target,intent(in) :: psps 460 !arrays 461 integer,intent(in) :: istwfk(nkpt),nband(nkpt*nsppol) 462 integer,intent(in) :: ngfft(18),nloalg(3),npwarr(nkpt),kg(3,mpw*mkmem),typat(natom) 463 integer,intent(in),target :: atindx(natom),nattyp(ntypat) 464 real(dp),intent(in) :: cg(2,mcg) 465 real(dp),intent(in) :: gmet(3,3),gprimd(3,3),kpt(3,nkpt),rmet(3,3) 466 real(dp),intent(in) :: xred(3,natom),ylm(:,:),ylmgr(:,:,:) 467 real(dp),intent(in),target :: ph1d(2,3*(2*mgfft+1)*natom) 468 type(pawcprj_type),intent(inout) :: cprj(ncprj,mcprj) 469 470 !Local variables------------------------------- 471 !scalars 472 integer :: blocksz,cg_bandpp,counter,cpopt,cprj_bandpp,dimffnl,ia,iatm,iatom1,iatom2 473 integer :: iband_max,iband_min,iband_start,ibg,ibgb,iblockbd,ibp,icg,icgb,icp1,icp2 474 integer :: ider,idir0,iend,ierr,ig,ii,ikg,ikpt,ilm,ipw,isize,isppol,istart,istwf_k,itypat,iwf1,iwf2,jdir 475 integer :: matblk,me_distrb,my_nspinor,n1,n1_2p1,n2,n2_2p1,n3,n3_2p1,kk,nlmn 476 integer :: mband,mband_cg,mband_cprj,npband_dfpt 477 integer :: nband_k,nband_cprj_k,nblockbd,ncpgr,nkpg,npband_bandfft,npws,npw_k,npw_nk,ntypat0 478 integer :: nband_cg_k 479 integer :: shift1,shift1b,shift2,shift2b,shift3,shift3b 480 integer :: spaceComm,spaceComm_band,spaceComm_fft,useylmgr 481 logical :: cg_band_distributed,cprj_band_distributed,one_atom 482 real(dp) :: arg 483 character(len=500) :: msg 484 !arrays 485 integer,allocatable :: bufsize(:),bufsize_wf(:),bufdisp(:),bufdisp_wf(:) 486 integer,allocatable :: dimlmn(:),kg_k(:,:),kg_k_loc(:,:) 487 integer,allocatable :: npw_block(:),npw_disp(:) 488 integer,pointer :: atindx_atm(:),indlmn_atm(:,:,:),nattyp_atm(:),pspso_atm(:) 489 real(dp) :: kpoint(3),work(6),tsec(2) 490 real(dp),allocatable :: cwavef(:,:),cwavef_tmp(:,:) 491 real(dp),allocatable :: ffnl(:,:,:,:),ffnl_npw(:,:,:,:),ffnl_tmp(:,:,:,:),ffnl_tmp_npw(:,:,:,:) 492 real(dp),allocatable :: kpg_k(:,:) 493 real(dp),allocatable :: ph3d(:,:,:),ph3d_npw(:,:,:),ph3d_tmp(:,:,:),ph3d_tmp_npw(:,:,:) 494 real(dp),allocatable :: phkxred(:,:),ylm_k(:,:),ylmgr_k(:,:,:) 495 real(dp),ABI_CONTIGUOUS pointer :: ekb_atm(:,:),ffspl_atm(:,:,:,:),ph1d_atm(:,:) 496 type(pawcprj_type),allocatable :: cwaveprj(:,:) 497 498 ! ********************************************************************* 499 500 DBG_ENTER('COLL') 501 502 !Nothing to do if current MPI process does treat kpoints or plane-waves 503 if (mcg==0.or.mcprj==0) return 504 505 !Preliminary tests 506 if (psps%useylm==0) then 507 msg='Not available for useylm=0!' 508 ABI_ERROR(msg) 509 end if 510 if ((choice<1.or.choice>6).and.choice/=23.and.choice/=24) then 511 msg='Bad choice!' 512 ABI_BUG(msg) 513 end if 514 if (idir>0.and.choice/=2.and.choice/=3.and.choice/=5) then 515 msg='Does not support idir>0 for that choice!' 516 ABI_BUG(msg) 517 end if 518 if (size(ylm)/=mpw*mkmem*mpsang*mpsang) then 519 msg='Wrong size for Ylm!' 520 ABI_BUG(msg) 521 end if 522 useylmgr=merge(0,1,(size(ylmgr)==0)) 523 if (useylmgr==0.and.(choice==3.or.choice==5.or.choice==23)) then 524 msg=' Ylm gradients have to be in memory for choice=3, 5, or 23!' 525 ABI_BUG(msg) 526 end if 527 528 !Init parallelism 529 npband_dfpt = 1 530 if (paral_kgb==1) then 531 me_distrb=mpi_enreg%me_kpt 532 spaceComm=mpi_enreg%comm_kpt 533 spaceComm_fft=mpi_enreg%comm_fft 534 npband_bandfft=mpi_enreg%nproc_band 535 cg_bandpp=mpi_enreg%bandpp 536 cprj_bandpp=mpi_enreg%bandpp 537 spaceComm_band=mpi_enreg%comm_band 538 cg_band_distributed=.true. 539 cprj_band_distributed=(mpi_enreg%nproc_band>1) 540 npband_dfpt=1 541 else 542 me_distrb=mpi_enreg%me_kpt 543 spaceComm=mpi_enreg%comm_cell 544 spaceComm_fft=xmpi_comm_self 545 npband_bandfft=1 546 cg_bandpp=1 547 cprj_bandpp=1 548 if (mpi_enreg%paralbd==1) then 549 spaceComm_band=mpi_enreg%comm_band ! not actually used as npband_bandfft=1 550 cg_band_distributed=.true. 551 cprj_band_distributed=.true. 552 !npband_dfpt=1 553 npband_dfpt=mpi_enreg%nproc_band 554 else 555 spaceComm_band=xmpi_comm_self 556 cg_band_distributed=.false. 557 cprj_band_distributed=.false. 558 end if 559 end if 560 if (cg_bandpp/=cprj_bandpp) then 561 ABI_BUG('cg_bandpp must be equal to cprj_bandpp!') 562 end if 563 564 !Manage parallelization over bands 565 mband=maxval(nband(1:nkpt*my_nsppol)) 566 mband_cg=mband/npband_bandfft/npband_dfpt 567 mband_cprj=mband_cg 568 569 !Manage parallelization over spinors 570 my_nspinor=max(1,nspinor/mpi_enreg%nproc_spinor) 571 572 !Check sizes 573 if (mpw*mband_cg*my_nspinor*mkmem*my_nsppol>mcg) then 574 ABI_BUG(sjoin('Bad mcg value!', itoa(mcg))) 575 end if 576 if (mband_cprj*my_nspinor*mkmem*my_nsppol>mcprj) then 577 ABI_BUG(sjoin('Bad mcprj value!', itoa(mcprj))) 578 end if 579 580 !Check sizes for cprj (distribution is tricky) 581 one_atom=(iatom>0) 582 if (one_atom.and.ncprj/=1) then 583 ABI_BUG('Bad value for ncprj dimension (should be 1) !') 584 end if 585 if (.not.one_atom.and.ncprj/=natom) then 586 ABI_BUG('Bad value for ncprj dimension (should be natom) !') 587 end if 588 589 !Initialize some variables 590 n1=ngfft(1);n2=ngfft(2);n3=ngfft(3) 591 n1_2p1=2*n1+1;n2_2p1=2*n2+1;n3_2p1=2*n3+1 592 ibg=0;icg=0;cpopt=0 593 ider=0;idir0=0;istart=idir;iend=idir 594 if (choice==3.or.choice==5.or.choice==23) ider=1 595 if (idir>0) then 596 if (choice==3) idir0=-idir 597 if (choice==5) idir0=idir 598 else 599 ! if (choice==23) idir0=-7 600 if (choice==3) idir0=-7 601 if (choice==5) idir0=4 602 end if 603 if (idir0==0.or.idir0==4) then 604 dimffnl=1+3*ider 605 else if (idir0/=-7) then 606 dimffnl=1+ider 607 else 608 dimffnl=1+6*ider 609 if(choice==3)then 610 istart=ider;iend=6*ider 611 end if 612 end if 613 nkpg=0 614 if (choice==3.or.choice==2.or.choice==23) nkpg=3*nloalg(3) 615 if (choice==4.or.choice==24) nkpg=9*nloalg(3) 616 617 !Set number of gradients for <p_i|Cnk> 618 ncpgr=0 619 if (idir==0) then 620 if (choice==2) ncpgr=3 621 if (choice==3) ncpgr=6 622 if (choice==23)ncpgr=9 623 if (choice==4) ncpgr=6 624 if (choice==24)ncpgr=9 625 if (choice==5) ncpgr=3 626 if (choice==6) ncpgr=63 627 else 628 ncpgr=1 629 end if 630 !Test cprj gradients dimension (just to be sure) 631 if (cprj(1,1)%ncpgr/=ncpgr) then 632 ABI_BUG('cprj are badly allocated !') 633 end if 634 635 636 !Extract data for treated atom(s) 637 if (one_atom) then 638 iatom1=iatom;iatom2=iatom 639 ntypat0=1;itypat=typat(iatom) 640 ABI_MALLOC(nattyp_atm,(ntypat0)) 641 nattyp_atm(1)=1 642 ABI_MALLOC(atindx_atm,(ntypat0)) 643 atindx_atm(1)=atindx(iatom) 644 ABI_MALLOC(ph1d_atm,(2,(n1_2p1+n2_2p1+n3_2p1)*ntypat0)) 645 shift1=(atindx(iatom)-1)*n1_2p1 646 shift2=(atindx(iatom)-1)*n2_2p1+natom*n1_2p1 647 shift3=(atindx(iatom)-1)*n3_2p1+natom*(n1_2p1+n2_2p1) 648 shift1b=0;shift2b=n1_2p1;shift3b=n1_2p1+n2_2p1 649 ph1d_atm(:,shift1b+1:shift1b+n1_2p1)=ph1d(:,shift1+1:shift1+n1_2p1) 650 ph1d_atm(:,shift2b+1:shift2b+n2_2p1)=ph1d(:,shift2+1:shift2+n2_2p1) 651 ph1d_atm(:,shift3b+1:shift3b+n3_2p1)=ph1d(:,shift3+1:shift3+n3_2p1) 652 ABI_MALLOC(ekb_atm,(psps%dimekb,ntypat0)) 653 ABI_MALLOC(indlmn_atm,(6,psps%lmnmax,ntypat0)) 654 ABI_MALLOC(ffspl_atm,(psps%mqgrid_ff,2,psps%lnmax,ntypat0)) 655 ABI_MALLOC(pspso_atm,(ntypat0)) 656 ekb_atm(:,1)=psps%ekb(:,itypat) 657 indlmn_atm(:,:,1)=psps%indlmn(:,:,itypat) 658 ffspl_atm(:,:,:,1)=psps%ffspl(:,:,:,itypat) 659 pspso_atm(1)=psps%pspso(itypat) 660 else 661 iatom1=1;iatom2=natom 662 ntypat0=ntypat 663 atindx_atm => atindx 664 nattyp_atm => nattyp 665 ph1d_atm => ph1d 666 ekb_atm => psps%ekb 667 indlmn_atm => psps%indlmn 668 ffspl_atm => psps%ffspl 669 pspso_atm => psps%pspso 670 end if 671 672 !Dimensioning and allocation of <p_i|Cnk> 673 ABI_MALLOC(dimlmn,(ncprj)) 674 dimlmn=0 ! Type-sorted cprj 675 if (one_atom) then 676 itypat=typat(iatom) 677 dimlmn(ia+1:ia+nattyp(itypat))=count(indlmn_atm(3,:,itypat)>0) 678 else 679 ia=0 680 do itypat=1,ntypat0 681 dimlmn(ia+1:ia+nattyp(itypat))=count(indlmn_atm(3,:,itypat)>0) 682 ia=ia+nattyp(itypat) 683 end do 684 end if 685 ABI_MALLOC(cwaveprj,(ncprj,my_nspinor*cprj_bandpp)) 686 call pawcprj_alloc(cwaveprj,ncpgr,dimlmn) 687 688 !Additional statements if band-fft parallelism 689 if (npband_bandfft>1) then 690 ABI_MALLOC(npw_block,(npband_bandfft)) 691 ABI_MALLOC(npw_disp,(npband_bandfft)) 692 !FB ABI_MALLOC(bufsize,(npband_bandfft*cg_bandpp)) 693 ABI_MALLOC(bufsize,(npband_bandfft)) 694 !FB ABI_MALLOC(bufdisp,(npband_bandfft*cg_bandpp)) 695 ABI_MALLOC(bufdisp,(npband_bandfft)) 696 !FB ABI_MALLOC(bufsize_wf,(npband_bandfft*cg_bandpp)) 697 ABI_MALLOC(bufsize_wf,(npband_bandfft)) 698 !FB ABI_MALLOC(bufdisp_wf,(npband_bandfft*cg_bandpp)) 699 ABI_MALLOC(bufdisp_wf,(npband_bandfft)) 700 end if 701 702 !Set output datastructure to zero 703 call pawcprj_set_zero(cprj) 704 705 !LOOP OVER SPINS 706 do isppol=1,my_nsppol 707 ikg=0 708 709 ! BIG FAT k POINT LOOP 710 do ikpt=1,nkpt 711 counter=100*ikpt+isppol 712 713 ! Select k point to be treated by this proc 714 nband_k=nband(ikpt+(isppol-1)*nkpt) 715 if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,isppol,me_distrb)) cycle 716 717 ! Retrieve k-point 718 kpoint(:)=kpt(:,ikpt) 719 istwf_k=istwfk(ikpt) 720 721 ! Retrieve number of plane waves 722 npw_k=npwarr(ikpt) 723 if (npband_bandfft>1) then 724 ! Special treatment for band-fft // 725 call xmpi_allgather(npw_k,npw_block,spaceComm_band,ierr) 726 npw_nk=sum(npw_block);npw_disp(1)=0 727 do ii=2,npband_bandfft 728 npw_disp(ii)=npw_disp(ii-1)+npw_block(ii-1) 729 end do 730 else 731 npw_nk=npw_k 732 end if 733 734 ! Retrieve (k+G) points and spherical harmonics 735 ABI_MALLOC(ylm_k,(npw_k,mpsang*mpsang)) 736 ABI_MALLOC(ylmgr_k,(npw_k,3,mpsang*mpsang*useylmgr)) 737 ABI_MALLOC(kg_k,(3,npw_nk)) 738 if (npband_bandfft>1) then 739 ! Special treatment for band-fft // 740 ABI_MALLOC(kg_k_loc,(3,npw_k)) 741 kg_k_loc(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg) 742 bufsize(:)=3*npw_block(:);bufdisp(:)=3*npw_disp(:) 743 call xmpi_allgatherv(kg_k_loc,3*npw_k,kg_k,bufsize,bufdisp,spaceComm_band,ierr) 744 else 745 kg_k(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg) 746 end if 747 do ilm=1,mpsang*mpsang 748 ylm_k(1:npw_k,ilm)=ylm(1+ikg:npw_k+ikg,ilm) 749 if (useylmgr>0) ylmgr_k(1:npw_k,1:3,ilm)=ylmgr(1+ikg:npw_k+ikg,1:3,ilm) 750 end do 751 752 ! Compute (k+G) vectors 753 ABI_MALLOC(kpg_k,(npw_nk,nkpg)) 754 if (nkpg>0) then 755 call mkkpg(kg_k,kpg_k,kpoint,nkpg,npw_nk) 756 end if 757 ! Allocate and compute the arrays phkxred and ph3d 758 ABI_MALLOC(phkxred,(2,ncprj)) 759 do ia=iatom1,iatom2 760 iatm=min(atindx_atm(ia),ncprj) 761 arg=two_pi*(kpoint(1)*xred(1,ia)+kpoint(2)*xred(2,ia)+kpoint(3)*xred(3,ia)) 762 phkxred(1,iatm)=cos(arg);phkxred(2,iatm)=sin(arg) 763 end do 764 matblk=ncprj;if (nloalg(2)<=0) matblk=0 765 ABI_MALLOC(ph3d,(2,npw_nk,matblk)) 766 if (matblk>0)then 767 ! Here, precomputation of ph3d 768 if (npband_bandfft>1) then 769 ! Special treatment for band-fft // 770 ABI_MALLOC(ph3d_tmp,(2,npw_k,matblk)) 771 call ph1d3d(1,ncprj,kg_k_loc,matblk,ncprj,npw_k,n1,n2,n3,phkxred,ph1d_atm,ph3d_tmp) 772 ABI_MALLOC(ph3d_tmp_npw,(2,matblk,npw_k)) 773 ABI_MALLOC(ph3d_npw,(2,matblk,npw_nk)) 774 isize=2*matblk;bufsize(:)=isize*npw_block(:);bufdisp(:)=isize*npw_disp(:) 775 do ipw=1,npw_k 776 ph3d_tmp_npw(:,:,ipw)=ph3d_tmp(:,ipw,:) 777 end do 778 call xmpi_allgatherv(ph3d_tmp_npw,isize*npw_k,ph3d_npw,bufsize,bufdisp,spaceComm_band,ierr) 779 do ipw=1,npw_nk 780 ph3d(:,ipw,:)=ph3d_npw(:,:,ipw) 781 end do 782 ABI_FREE(ph3d_npw) 783 ABI_FREE(ph3d_tmp_npw) 784 ABI_FREE(ph3d_tmp) 785 else 786 call ph1d3d(1,ncprj,kg_k,matblk,ncprj,npw_k,n1,n2,n3,phkxred,ph1d_atm,ph3d) 787 end if 788 else if (npband_bandfft>1) then 789 ABI_ERROR('Band-fft parallelism +nloag(1)<0 forbidden !') 790 end if 791 792 ! Compute nonlocal form factors ffnl at all (k+G) 793 ABI_MALLOC(ffnl,(npw_nk,dimffnl,psps%lmnmax,ntypat0)) 794 if (npband_bandfft>1) then 795 ! Special treatment for band-fft // 796 ABI_MALLOC(ffnl_tmp,(npw_k,dimffnl,psps%lmnmax,ntypat0)) 797 call mkffnl(psps%dimekb,dimffnl,ekb_atm,ffnl_tmp,ffspl_atm,& 798 & gmet,gprimd,ider,idir0,indlmn_atm,kg_k_loc,kpg_k,kpoint,psps%lmnmax,& 799 & psps%lnmax,psps%mpsang,psps%mqgrid_ff,nkpg,npw_k,ntypat0,& 800 & pspso_atm,psps%qgrid_ff,rmet,psps%usepaw,psps%useylm,ylm_k,ylmgr_k) 801 ABI_MALLOC(ffnl_tmp_npw,(dimffnl,psps%lmnmax,ntypat0,npw_k)) 802 ABI_MALLOC(ffnl_npw,(dimffnl,psps%lmnmax,ntypat0,npw_nk)) 803 isize=dimffnl*psps%lmnmax*ntypat0 804 bufsize(:)=isize*npw_block(:);bufdisp(:)=isize*npw_disp(:) 805 do ipw=1,npw_k 806 ffnl_tmp_npw(:,:,:,ipw)=ffnl_tmp(ipw,:,:,:) 807 end do 808 call xmpi_allgatherv(ffnl_tmp_npw,isize*npw_k,ffnl_npw,bufsize,bufdisp,spaceComm_band,ierr) 809 do ipw=1,npw_nk 810 ffnl(ipw,:,:,:)=ffnl_npw(:,:,:,ipw) 811 end do 812 ABI_FREE(ffnl_npw) 813 ABI_FREE(ffnl_tmp_npw) 814 ABI_FREE(ffnl_tmp) 815 else 816 call mkffnl(psps%dimekb,dimffnl,ekb_atm,ffnl,ffspl_atm,& 817 & gmet,gprimd,ider,idir0,indlmn_atm,kg_k,kpg_k,kpoint,psps%lmnmax,& 818 & psps%lnmax,psps%mpsang,psps%mqgrid_ff,nkpg,npw_k,ntypat0,& 819 & pspso_atm,psps%qgrid_ff,rmet,psps%usepaw,psps%useylm,ylm_k,ylmgr_k) 820 end if 821 822 ! No more need of kg_g_tmp 823 if (npband_bandfft>1) then 824 ABI_FREE(kg_k_loc) 825 end if 826 827 ! Allocate arrays for a wave-function (or a block of WFs) 828 ABI_MALLOC(cwavef,(2,npw_nk*my_nspinor*cg_bandpp)) 829 if (npband_bandfft>1) then 830 isize=2*my_nspinor*cg_bandpp;bufsize(:)=isize*npw_block(:);bufdisp(:)=isize*npw_disp(:) 831 isize=2*my_nspinor*npw_k*cg_bandpp;bufsize_wf(:)=isize 832 !FB do ii=1,npband_bandfft*cg_bandpp 833 do ii=1,npband_bandfft 834 bufdisp_wf(ii)=(ii-1)*isize 835 end do 836 end if 837 838 ! Loop over bands or blocks of bands 839 icgb=icg ; ibgb=ibg ; iband_start=1 840 blocksz=npband_bandfft*cg_bandpp 841 nblockbd=nband_k/blocksz 842 nband_cprj_k=merge(nband_k/npband_bandfft/npband_dfpt,nband_k,cprj_band_distributed) 843 nband_cg_k=merge(nband_k/npband_bandfft/npband_dfpt,nband_k,cg_band_distributed) 844 do iblockbd=1,nblockbd 845 iband_min=1+(iblockbd-1)*blocksz 846 iband_max=iblockbd*blocksz 847 848 if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,iband_min,iband_max,isppol,me_distrb)) then 849 if (.not.cg_band_distributed) icgb=icgb+npw_k*my_nspinor*blocksz 850 if (.not.cprj_band_distributed) ibgb=ibgb+my_nspinor*blocksz 851 cycle 852 end if 853 854 ! Extract wavefunction information 855 ! Special treatment for band-fft parallelism 856 if (npband_bandfft>1) then 857 !Transpose WF to get them in "FFT" representation 858 ABI_MALLOC(cwavef_tmp,(2,npw_k*my_nspinor*blocksz)) 859 do ig=1,npw_k*my_nspinor*blocksz 860 cwavef_tmp(1,ig)=cg(1,ig+icgb) 861 cwavef_tmp(2,ig)=cg(2,ig+icgb) 862 end do 863 call xmpi_alltoallv(cwavef_tmp,bufsize_wf,bufdisp_wf,cwavef,bufsize,bufdisp,spaceComm_band,ierr) 864 ABI_FREE(cwavef_tmp) 865 !Reorder WF according to cg_bandpp and/or spinor 866 if (cg_bandpp>1.or.my_nspinor>1) then 867 !FB ABI_MALLOC(cwavef_tmp,(2,npw_nk*my_nspinor*blocksz)) 868 ABI_MALLOC(cwavef_tmp,(2,npw_nk*my_nspinor*cg_bandpp)) 869 !FB do ig=1,npw_nk*my_nspinor*blocksz 870 do ig=1,npw_nk*my_nspinor*cg_bandpp 871 cwavef_tmp(:,ig)=cwavef(:,ig) 872 end do 873 shift1=0 874 do iwf2=1,cg_bandpp 875 do ig=1,my_nspinor 876 shift2=0 877 do iwf1=1,npband_bandfft 878 npws=npw_block(iwf1) 879 ipw=shift2+(iwf2-1)*my_nspinor*npws+(ig-1)*npws 880 cwavef(:,shift1+1:shift1+npws)=cwavef_tmp(:,ipw+1:ipw+npws) 881 shift1=shift1+npws ; shift2=shift2+cg_bandpp*my_nspinor*npws 882 end do 883 end do 884 end do 885 ABI_FREE(cwavef_tmp) 886 end if 887 else 888 do ig=1,npw_k*my_nspinor*cg_bandpp 889 cwavef(1,ig)=cg(1,ig+icgb) 890 cwavef(2,ig)=cg(2,ig+icgb) 891 end do 892 end if 893 894 ! Compute scalar product of wavefunction with all NL projectors 895 do ibp=1,cg_bandpp ! Note: we suppose cp_bandpp=cprj_bandpp 896 iwf1=1+(ibp-1)*npw_nk*my_nspinor;iwf2=ibp*npw_nk*my_nspinor 897 icp1=1+(ibp-1)*my_nspinor;icp2=ibp*my_nspinor 898 do jdir=istart,iend 899 call timab(1294,1,tsec) 900 call getcprj(choice,cpopt,cwavef(:,iwf1:iwf2),cwaveprj(:,icp1:icp2),& 901 & ffnl,jdir,indlmn_atm,istwf_k,kg_k,kpg_k,kpoint,psps%lmnmax,& 902 & mgfft,mpi_enreg,ncprj,nattyp_atm,ngfft,nloalg,& 903 & npw_nk,my_nspinor,ntypat0,phkxred,ph1d_atm,ph3d,ucvol,psps%useylm) 904 call timab(1294,2,tsec) 905 end do 906 end do 907 ! Export cwaveprj to big array cprj 908 call pawcprj_put(atindx_atm,cwaveprj,cprj,ncprj,iband_start,ibgb,ikpt,iorder_cprj,isppol,& 909 & mband_cprj,mkmem,natom,cprj_bandpp,nband_cprj_k,dimlmn,my_nspinor,nsppol,uncp,& 910 & mpi_comm_band=spaceComm_band,to_be_gathered=(cg_band_distributed.and.(.not.cprj_band_distributed))) 911 912 iband_start=iband_start+merge(cg_bandpp,blocksz,cprj_band_distributed) 913 914 ! End loop over bands 915 icgb=icgb+npw_k*my_nspinor*blocksz 916 end do 917 918 ! Shift array memory (if mkmem/=0) 919 if (mkmem/=0) then 920 ibg=ibg+my_nspinor*nband_cprj_k 921 !FB icg=icg+my_nspinor*nband_cg_k*npw_k 922 icg=icg+my_nspinor*nband_k*npw_k 923 ikg=ikg+npw_k 924 end if 925 926 ! End big k point loop 927 ABI_FREE(ffnl) 928 ABI_FREE(ph3d) 929 ABI_FREE(phkxred) 930 ABI_FREE(kg_k) 931 ABI_FREE(kpg_k) 932 ABI_FREE(ylm_k) 933 ABI_FREE(ylmgr_k) 934 ABI_FREE(cwavef) 935 end do 936 ! End loop over spins 937 end do 938 939 if ((iatom<=0).and.(choice==23)) then 940 do iatom1=1,ncprj 941 do ii=1,mcprj 942 nlmn=cprj(iatom1,ii)%nlmn 943 do kk=1,nlmn 944 work(1:6)=cprj(iatom1,ii)%dcp(1,1:6,kk) 945 call strconv(work,gprimd,work) 946 cprj(iatom1,ii)%dcp(1,1:6,kk)=work(1:6) 947 work(1:6)=cprj(iatom1,ii)%dcp(2,1:6,kk) 948 call strconv(work,gprimd,work) 949 cprj(iatom1,ii)%dcp(2,1:6,kk)=work(1:6) 950 end do 951 end do 952 end do 953 end if 954 955 !If needed, gather computed scalars 956 if (.not.(cg_band_distributed .and. cprj_band_distributed)) then 957 call pawcprj_mpi_sum(cprj,spaceComm_band,ierr) 958 end if 959 960 !Deallocate temporary storage 961 if (one_atom) then 962 ABI_FREE(atindx_atm) 963 ABI_FREE(nattyp_atm) 964 ABI_FREE(ph1d_atm) 965 ABI_FREE(ekb_atm) 966 ABI_FREE(indlmn_atm) 967 ABI_FREE(ffspl_atm) 968 ABI_FREE(pspso_atm) 969 end if 970 nullify(atindx_atm,nattyp_atm,ph1d_atm,ekb_atm,indlmn_atm,ffspl_atm,pspso_atm) 971 call pawcprj_free(cwaveprj) 972 ABI_FREE(cwaveprj) 973 ABI_FREE(dimlmn) 974 if (npband_bandfft>1) then 975 ABI_FREE(npw_block) 976 ABI_FREE(npw_disp) 977 ABI_FREE(bufsize) 978 ABI_FREE(bufdisp) 979 ABI_FREE(bufsize_wf) 980 ABI_FREE(bufdisp_wf) 981 end if 982 983 DBG_EXIT('COLL') 984 985 end subroutine ctocprj
ABINIT/getcprj [ Functions ]
NAME
getcprj
FUNCTION
Compute <Proj_i|Cnk> for one wave function |Cnk> expressed in reciprocal space. Compute also derivatives of <Proj_i|Cnk>. |Proj_i> are non-local projectors (for each atom and each l,m,n)
INPUTS
choice=chooses possible output: In addition to projected wave function: choice=1 => nothing else =2 => 1st gradients with respect to atomic position(s) =3 => 1st gradients with respect to strain(s) =23=> 1st gradients with respect to strain(s) and atm pos =4 => 2nd derivatives with respect to atomic pos. =24=> 1st and 2nd derivatives with respect to atomic pos. =5 => 1st gradients with respect to k wavevector =6 => 2nd derivatives with respect to strain and atm. pos. cpopt=1 if <Proj_i|Cnk> are already in memory; see below (side effects). cwavef(2,nspinor*npw_k)=input cmplx wavefunction coefficients <G|Cnk> ffnl(npw_k,dimffnl,lmnmax,ntypat)=nonlocal form factors to be used for the application of the nl operator idir=direction of the derivative, i.e. dir. of - atom to be moved in the case choice=2 - strain component in the case choice=3 - k point direction in the case choice=5 Compatible only with choice=2,3,5; if idir=0, all derivatives are computed indlmn(6,i,ntypat)= array giving l,m,n,lm,ln,s for i=lmn istwf_k=option parameter that describes the storage of wfs kg_k(3,npw_k)=reduced planewave coordinates kpg(npw_k,npkg)=(k+G) components and related data kpoint(3)=k point in terms of recip. translations lmnmax=max. number of (l,m,n) components over all types of atoms mgfft=maximum size of 1D FFTs mpi_enreg=information about MPI parallelization natom=number of atoms in cell nattyp(ntypat)=number of atoms of each type ngfft(18)=contain all needed information about 3D FFT, see ~ABINIT/Infos/vargs.htm#ngfft nloalg(3)=governs the choice of the algorithm for nonlocal operator npw_k=number of planewaves for given k point nspinor=number of spinorial components of the wavefunctions (on current proc) ntypat=number of types of atoms in unit cell phkxred(2,natom)=phase factors exp(2 pi kpoint.xred) ph1d(2,3*(2*mgfft+1)*natom)=1-dim structure factor phase information ph3d(2,npw_k,natom)=3D structure factors, for each atom and plane wave only used if nloalg(2)>0 ucvol= unit cell volume useylm=governs the way the nonlocal operator is to be applied
SIDE EFFECTS
cwaveprj(natom,nspinor) <type(pawcprj_type)>=projected input wave function <Proj_i|Cnk> with all NL projectors (and derivatives) if cpopt=1 the projected scalars have been already been computed and only derivatives are computed here if cpopt=0 the projected scalars and derivatives are computed here
TODO
Spin-orbit
SOURCE
117 subroutine getcprj(choice,cpopt,cwavef,cwaveprj,ffnl,& 118 & idir,indlmn,istwf_k,kg_k,kpg,kpoint,lmnmax,mgfft,mpi_enreg,& 119 & natom,nattyp,ngfft,nloalg,npw_k,nspinor,ntypat,& 120 & phkxred,ph1d,ph3d,ucvol,useylm) 121 122 !Arguments ------------------------------- 123 !scalars 124 integer,intent(in) :: choice,cpopt,idir,istwf_k,lmnmax 125 integer,intent(in) :: mgfft,natom,npw_k,nspinor,ntypat,useylm 126 real(dp),intent(in) :: ucvol 127 type(MPI_type),intent(in) :: mpi_enreg 128 !arrays 129 integer,intent(in) :: indlmn(6,lmnmax,ntypat),kg_k(3,npw_k),nattyp(ntypat) 130 integer,intent(in) :: ngfft(18),nloalg(3) 131 real(dp),intent(in) :: cwavef(2,npw_k*nspinor) 132 real(dp),intent(in),target :: ffnl(:,:,:,:),kpg(:,:),ph3d(:,:,:) 133 real(dp),intent(in) :: kpoint(3),ph1d(2,3*(2*mgfft+1)*natom),phkxred(2,natom) 134 type(pawcprj_type),intent(inout) :: cwaveprj(natom,nspinor) 135 136 !Local variables------------------------------- 137 !scalars 138 logical :: no_opernla_mv 139 integer :: choice_,cplex,dimffnl,ia,ia1,ia2,ia3,ia4,iatm,ic,ii,ilmn,ishift,ispinor,itypat 140 integer :: jc,matblk,mincat,nd2gxdt,ndgxdt,nincat,nkpg,nkpg_,nlmn,signs 141 !arrays 142 real(dp) :: tsec(2) 143 integer,allocatable :: cplex_dgxdt(:),cplex_d2gxdt(:),indlmn_typ(:,:) 144 real(dp),allocatable :: d2gxdt(:,:,:,:,:),dgxdt(:,:,:,:,:),ffnl_typ(:,:,:) 145 real(dp),allocatable :: gx(:,:,:,:) 146 real(dp), pointer :: kpg_(:,:),ph3d_(:,:,:) 147 148 ! ********************************************************************* 149 150 DBG_ENTER('COLL') 151 152 call timab(1290,1,tsec) 153 154 !Nothing to do in that case 155 if (cpopt==1.and.choice==1) return 156 157 !Not available for useylm=0 158 if (useylm==0) then 159 ABI_ERROR('Not available for useylm=0 !') 160 end if 161 162 !Error on bad choice 163 if ((choice<1.or.choice>6).and.choice/=23.and.choice/=24) then 164 ABI_BUG('Does not presently support this choice !') 165 end if 166 167 !Error on bad idir 168 if (idir>0.and.choice/=2.and.choice/=3.and.choice/=5) then 169 ABI_BUG('Does not support idir>0 for this choice') 170 end if 171 172 !Error on sizes 173 nkpg=size(kpg,2) 174 if (nkpg>0) then 175 if( (choice==2.and.nkpg<3) .or. & 176 & ((choice==4.or.choice==24).and.nkpg<9) .or. & 177 & ((choice==6.or.choice==3.or.choice==23).and.nkpg<3) ) then 178 ABI_BUG('Incorrect size for kpg array !') 179 end if 180 end if 181 if (size(ffnl,1)/=npw_k.or.size(ffnl,3)/=lmnmax) then 182 ABI_BUG('Incorrect size for ffnl!') 183 end if 184 if (size(ph3d)>0) then 185 if (size(ph3d,2)/=npw_k) then 186 ABI_BUG('Incorrect size for ph3d!') 187 end if 188 end if 189 190 no_opernla_mv = nloalg(1)==4.or.nloalg(1)==8.or.nloalg(1)==10 ! have to be consistent with nonlop_ylm 191 192 !Define dimensions of projected scalars 193 dimffnl=size(ffnl,2) 194 ndgxdt=0;nd2gxdt=0 195 if (idir==0) then 196 if (choice==2) ndgxdt=3 197 if (choice==3) ndgxdt=6 198 if (choice==23) ndgxdt=9 199 if (choice==4) nd2gxdt=6 200 if (choice==24) then 201 ndgxdt=3;nd2gxdt=6 202 end if 203 if (choice==5) ndgxdt=3 204 if (choice==6) then 205 ndgxdt=9;nd2gxdt=54 206 end if 207 else 208 ndgxdt=1 209 end if 210 if(cwaveprj(1,1)%ncpgr<ndgxdt+nd2gxdt) then 211 ABI_BUG('Incorrect size for ncpgr') 212 end if 213 !Eventually re-compute (k+G) vectors (and related data) 214 if (nkpg==0) then 215 nkpg_=0 216 if (choice==4.or.choice==24) nkpg_=9 217 if (choice==2.or.choice==3.or.choice==23) nkpg_=3 218 ABI_MALLOC(kpg_,(npw_k,nkpg_)) 219 if (nkpg_>0) then 220 call mkkpg(kg_k,kpg_,kpoint,nkpg_,npw_k) 221 end if 222 else 223 nkpg_=nkpg 224 kpg_ => kpg 225 end if 226 227 !Some other dims 228 mincat=min(NLO_MINCAT,maxval(nattyp)) 229 cplex=2;if (istwf_k>1) cplex=1 230 choice_=choice;if (cpopt==1) choice_=-choice 231 signs=1;if (idir>0) signs=2 232 !Eventually allocate temporary array for ph3d 233 if (nloalg(2)<=0) then 234 matblk=mincat 235 ABI_MALLOC(ph3d_,(2,npw_k,matblk)) 236 else 237 matblk=size(ph3d,3) 238 ph3d_ => ph3d 239 end if 240 241 !Loop over atom types 242 ia1=1;iatm=0 243 do itypat=1,ntypat 244 ia2=ia1+nattyp(itypat)-1;if (ia2<ia1) cycle 245 nlmn=count(indlmn(3,:,itypat)>0) 246 247 ! Retrieve some data for this type of atom 248 ABI_MALLOC(indlmn_typ,(6,nlmn)) 249 ABI_MALLOC(ffnl_typ,(npw_k,dimffnl,nlmn)) 250 indlmn_typ(:,1:nlmn)=indlmn(:,1:nlmn,itypat) 251 ffnl_typ(:,:,1:nlmn)=ffnl(:,:,1:nlmn,itypat) 252 253 ! Loop on blocks of atoms inside type 254 do ia3=ia1,ia2,mincat 255 ia4=min(ia2,ia3+mincat-1);nincat=ia4-ia3+1 256 ! Prepare the phase factors if they were not already computed 257 if (nloalg(2)<=0) then 258 call ph1d3d(ia3,ia4,kg_k,matblk,natom,npw_k,ngfft(1),ngfft(2),ngfft(3),& 259 & phkxred,ph1d,ph3d_) 260 end if 261 262 ! Allocate memory for projected scalars 263 ABI_MALLOC(gx,(cplex,nlmn,nincat,nspinor)) 264 ABI_MALLOC(dgxdt,(cplex,ndgxdt,nlmn,nincat,nspinor)) 265 ABI_MALLOC(d2gxdt,(cplex,nd2gxdt,nlmn,nincat,nspinor)) 266 ABI_MALLOC(cplex_dgxdt,(ndgxdt)) 267 ABI_MALLOC(cplex_d2gxdt,(nd2gxdt)) 268 269 ! Retrieve eventually <p_i|c> coeffs 270 if (cpopt==1) then 271 do ispinor=1,nspinor 272 do ia=1,nincat 273 gx(1:cplex,1:nlmn,ia,ispinor)=cwaveprj(iatm+ia,ispinor)%cp(1:cplex,1:nlmn) 274 end do 275 end do 276 end if 277 278 ! Compute <p_i|c> scalars (and derivatives) for this block of atoms 279 if (abs(choice_)>1.or.no_opernla_mv) then 280 call timab(1291,1,tsec) 281 call opernla_ylm(choice_,cplex,cplex_dgxdt,cplex_d2gxdt,dimffnl,d2gxdt,dgxdt,ffnl_typ,gx,& 282 & ia3,idir,indlmn_typ,istwf_k,kpg_,matblk,mpi_enreg,nd2gxdt,ndgxdt,nincat,nkpg_,nlmn,& 283 & nloalg,npw_k,nspinor,ph3d_,signs,ucvol,cwavef) 284 call timab(1291,2,tsec) 285 else 286 call timab(1292,1,tsec) 287 call opernla_ylm_mv(choice_,cplex,dimffnl,ffnl_typ,gx,& 288 & ia3,indlmn_typ,istwf_k,matblk,mpi_enreg,nincat,nlmn,& 289 & nloalg,npw_k,nspinor,ph3d_,ucvol,cwavef) 290 call timab(1292,2,tsec) 291 end if 292 293 ! Transfer result to output variable cwaveprj 294 if (cpopt==0) then 295 do ispinor=1,nspinor 296 do ia=1,nincat 297 cwaveprj(iatm+ia,ispinor)%nlmn=nlmn 298 cwaveprj(iatm+ia,ispinor)%cp(1:cplex,1:nlmn)=gx(1:cplex,1:nlmn,ia,ispinor) 299 if(cplex==1) cwaveprj(iatm+ia,ispinor)%cp(2,1:nlmn)=zero 300 end do 301 end do 302 end if 303 if (cpopt>=0.and.choice>1) then 304 ishift=0 305 if ((idir>0).and.(cwaveprj(1,1)%ncpgr>ndgxdt)) ishift=idir-1 306 if(cplex==2)then 307 do ispinor=1,nspinor 308 do ia=1,nincat 309 ! cwaveprj(iatm+ia,ispinor)%ncpgr=ndgxdt+nd2gxdt 310 if (ndgxdt>0) cwaveprj(iatm+ia,ispinor)%dcp(1:2,1+ishift:ndgxdt+ishift,1:nlmn)=& 311 & dgxdt(1:2,1:ndgxdt,1:nlmn,ia,ispinor) 312 if (nd2gxdt>0)cwaveprj(iatm+ia,ispinor)%dcp(1:2,ndgxdt+1+ishift:ndgxdt+nd2gxdt+ishift,1:nlmn)=& 313 & d2gxdt(1:2,1:nd2gxdt,1:nlmn,ia,ispinor) 314 end do 315 end do 316 else 317 ! cplex_dgxdt(i) = 1 if dgxdt(1,i,:,:) is real, 2 if it is pure imaginary 318 ! cplex_d2gxdt(i) = 1 if d2gxdt(1,i,:,:) is real, 2 if it is pure imaginary 319 do ispinor=1,nspinor 320 do ia=1,nincat 321 ! cwaveprj(iatm+ia,ispinor)%ncpgr=ndgxdt+nd2gxdt 322 if (ndgxdt>0) then 323 do ilmn =1,nlmn 324 do ii = 1,ndgxdt 325 ic = cplex_dgxdt(ii) ; jc = 3 - ic 326 cwaveprj(iatm+ia,ispinor)%dcp(ic,ii+ishift,ilmn)=dgxdt(1,ii,ilmn,ia,ispinor) 327 cwaveprj(iatm+ia,ispinor)%dcp(jc,ii+ishift,ilmn)=zero 328 end do 329 end do 330 end if 331 if (nd2gxdt>0) then 332 do ilmn =1,nlmn 333 do ii = 1,nd2gxdt 334 ic = cplex_d2gxdt(ii) ; jc = 3 - ic 335 cwaveprj(iatm+ia,ispinor)%dcp(ic,ndgxdt+ii+ishift,ilmn)=d2gxdt(1,ii,ilmn,ia,ispinor) 336 cwaveprj(iatm+ia,ispinor)%dcp(jc,ndgxdt+ii+ishift,ilmn)=zero 337 end do 338 end do 339 end if 340 end do 341 end do 342 end if 343 end if 344 345 ! End loop inside block of atoms 346 iatm=iatm+nincat 347 ABI_FREE(gx) 348 ABI_FREE(dgxdt) 349 ABI_FREE(d2gxdt) 350 ABI_FREE(cplex_dgxdt) 351 ABI_FREE(cplex_d2gxdt) 352 end do 353 354 ! End loop over atom types 355 ia1=ia2+1 356 ABI_FREE(indlmn_typ) 357 ABI_FREE(ffnl_typ) 358 end do 359 360 if (nkpg==0) then 361 ABI_FREE(kpg_) 362 end if 363 if (nloalg(2)<=0) then 364 ABI_FREE(ph3d_) 365 end if 366 367 call timab(1290,2,tsec) 368 369 DBG_EXIT('COLL') 370 371 end subroutine getcprj
ABINIT/m_cgprj [ Modules ]
NAME
m_cgprj
FUNCTION
Routines to compute <Proj_i|Cnk> with |Cnk> expressed in reciprocal space.
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (MT) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_cgprj 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 use m_xmpi 28 29 use defs_abitypes, only : MPI_type 30 use defs_datatypes, only : pseudopotential_type 31 use m_fstrings, only : itoa, sjoin 32 use m_kg, only : ph1d3d, mkkpg 33 use m_geometry, only : strconv 34 use m_mkffnl, only : mkffnl 35 use m_mpinfo, only : proc_distrb_cycle 36 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_put, pawcprj_free, & 37 pawcprj_set_zero, pawcprj_mpi_sum, pawcprj_copy, pawcprj_lincom 38 use m_opernla_ylm, only : opernla_ylm 39 use m_opernla_ylm_mv, only : opernla_ylm_mv 40 use m_time, only : timab 41 use m_io_tools, only : flush_unit 42 43 implicit none 44 45 private