TABLE OF CONTENTS
- ABINIT/bandfft_kpt_mpi_recv
- ABINIT/bandfft_kpt_mpi_send
- ABINIT/m_bandfft_kpt
- ABINIT/prep_bandfft_tabs
- m_bandfft_kpt/bandfft_kpt_copy
- m_bandfft_kpt/bandfft_kpt_destroy
- m_bandfft_kpt/bandfft_kpt_destroy_array
- m_bandfft_kpt/bandfft_kpt_get_ikpt
- m_bandfft_kpt/bandfft_kpt_init1
- m_bandfft_kpt/bandfft_kpt_init2
- m_bandfft_kpt/bandfft_kpt_reset
- m_bandfft_kpt/bandfft_kpt_restoretabs
- m_bandfft_kpt/bandfft_kpt_savetabs
- m_bandfft_kpt/bandfft_kpt_set_ikpt
- m_bandfft_kpt/bandfft_kpt_type
ABINIT/bandfft_kpt_mpi_recv [ Functions ]
NAME
bandfft_kpt_mpi_recv
FUNCTION
Receive a bandfft_kpt_type inside a MPI communicator.
INPUTS
sender=ID of the sender in spaceComm tag=message tag spaceComm=MPI Communicator
OUTPUT
ierr=Error status output=# of on proc. sender
SOURCE
1546 subroutine bandfft_kpt_mpi_recv(output,sender,tag,spaceComm,ierr) 1547 1548 !Arguments ------------------------------------ 1549 !scalars 1550 integer,intent(in) :: sender,spaceComm,tag 1551 integer,intent(out) :: ierr 1552 !arrays 1553 type(bandfft_kpt_type),intent(out) :: output 1554 1555 !Local variables------------------------------- 1556 !scalars 1557 integer :: ipck,nsize,size_dp,size_int 1558 integer :: size1_kg_k_gather=0,size2_kg_k_gather=0,size_recvcounts=0 1559 integer :: size_sendcounts=0,size_rdispls=0,size_sdispls=0,size1_gbound=0,size2_gbound 1560 integer :: size1_kg_k_gather_sym=0,size2_kg_k_gather_sym=0,size_rdispls_sym=0 1561 integer :: size_sdispls_sym=0,size_recvcounts_sym=0,size_recvcounts_sym_tot=0 1562 integer :: size_sendcounts_sym=0,size_sendcounts_sym_all=0,size_tab_proc=0 1563 integer :: size_indices_pw_fft=0,size_sendcount_fft=0,size_senddisp_fft=0 1564 integer :: size_recvcount_fft=0,size_recvdisp_fft=0,size1_kg_k_fft=0,size2_kg_k_fft=0 1565 integer :: size1_ffnl_gather=0,size2_ffnl_gather=0,size3_ffnl_gather=0,size4_ffnl_gather=0 1566 integer :: size_kinpw_gather=0,size1_ph3d_gather=0,size2_ph3d_gather=0,size3_ph3d_gather=0 1567 integer :: size1_kpg_k_gather=0,size2_kpg_k_gather=0,sz1,sz2,sz3,sz4 1568 !arrays 1569 integer,allocatable :: buffer_int(:) 1570 real(dp),allocatable :: buffer_dp(:) 1571 1572 ! ************************************************************************* 1573 1574 if (xmpi_comm_size(spaceComm)<=1) return 1575 1576 !=== Receive flags and array sizes ==== 1577 nsize=45 1578 ABI_MALLOC(buffer_int,(nsize)) 1579 call xmpi_recv(buffer_int,sender,3*tag-2,spaceComm,ierr) 1580 output%flag1_is_allocated=buffer_int(1) 1581 output%flag2_is_allocated=buffer_int(2) 1582 output%flag3_is_allocated=buffer_int(3) 1583 output%have_to_reequilibrate=(buffer_int(4)/=0) 1584 output%istwf_k=buffer_int(5) 1585 output%npw_tot=buffer_int(6) 1586 output%npw_fft=buffer_int(7) 1587 output%ndatarecv=buffer_int(8) 1588 output%idatarecv0=buffer_int(9) 1589 output%ndatarecv_tot=buffer_int(10) 1590 output%ndatasend_sym=buffer_int(11) 1591 size_recvcounts=buffer_int(12) 1592 size_sendcounts=buffer_int(13) 1593 size_rdispls=buffer_int(14) 1594 size_sdispls=buffer_int(15) 1595 size1_gbound=buffer_int(16) 1596 size2_gbound=buffer_int(17) 1597 size1_ffnl_gather=buffer_int(18) 1598 size2_ffnl_gather=buffer_int(19) 1599 size3_ffnl_gather=buffer_int(20) 1600 size4_ffnl_gather=buffer_int(21) 1601 size_kinpw_gather=buffer_int(22) 1602 size1_ph3d_gather=buffer_int(23) 1603 size2_ph3d_gather=buffer_int(24) 1604 size3_ph3d_gather=buffer_int(25) 1605 size1_kpg_k_gather=buffer_int(26) 1606 size2_kpg_k_gather=buffer_int(27) 1607 size_indices_pw_fft=buffer_int(28) 1608 size_sendcount_fft=buffer_int(29) 1609 size_senddisp_fft=buffer_int(30) 1610 size_recvcount_fft=buffer_int(31) 1611 size_recvdisp_fft=buffer_int(32) 1612 size1_kg_k_fft=buffer_int(33) 1613 size2_kg_k_fft=buffer_int(34) 1614 size1_kg_k_gather=buffer_int(35) 1615 size2_kg_k_gather=buffer_int(36) 1616 size1_kg_k_gather_sym=buffer_int(37) 1617 size2_kg_k_gather_sym=buffer_int(38) 1618 size_rdispls_sym=buffer_int(39) 1619 size_sdispls_sym=buffer_int(40) 1620 size_recvcounts_sym=buffer_int(41) 1621 size_recvcounts_sym_tot=buffer_int(42) 1622 size_sendcounts_sym=buffer_int(43) 1623 size_sendcounts_sym_all=buffer_int(44) 1624 size_tab_proc=buffer_int(45) 1625 ABI_FREE(buffer_int) 1626 1627 !=== Compute amount of transmitted data ==== 1628 size_int=size1_kg_k_gather*size2_kg_k_gather & 1629 & +size_recvcounts+size_sendcounts+size_rdispls+size_sdispls & 1630 & +size1_gbound*size2_gbound & 1631 & +size1_kg_k_gather_sym*size2_kg_k_gather_sym & 1632 & +size_rdispls_sym+size_sdispls_sym & 1633 & +size_recvcounts_sym+size_recvcounts_sym_tot & 1634 & +size_sendcounts_sym+size_sendcounts_sym_all+size_tab_proc & 1635 & +size_indices_pw_fft+size_sendcount_fft+size_senddisp_fft & 1636 & +size_recvcount_fft+size_recvdisp_fft & 1637 & +size1_kg_k_fft*size2_kg_k_fft 1638 size_dp=size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather & 1639 & +size_kinpw_gather & 1640 & +size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather & 1641 & +size1_kpg_k_gather*size2_kpg_k_gather 1642 1643 !=== Receive and unpack integer data ==== 1644 if (size_int>0) then 1645 ipck=0 1646 ABI_MALLOC(buffer_int,(size_int)) 1647 call xmpi_recv(buffer_int,sender,3*tag-1,spaceComm,ierr) 1648 1649 #if defined HAVE_GPU && defined HAVE_YAKL 1650 if (associated(output%kg_k_gather)) then 1651 if(output%gpu_option==ABI_GPU_KOKKOS) then 1652 ABI_FREE_MANAGED(output%kg_k_gather) 1653 else 1654 ABI_FREE(output%kg_k_gather) 1655 end if 1656 end if 1657 if (size1_kg_k_gather*size2_kg_k_gather>0) then 1658 nsize=size1_kg_k_gather*size2_kg_k_gather 1659 sz1=size1_kg_k_gather;sz2=size2_kg_k_gather 1660 if(output%gpu_option==ABI_GPU_KOKKOS) then 1661 ABI_MALLOC_MANAGED(output%kg_k_gather,(/sz1,sz2/)) 1662 else 1663 ABI_MALLOC(output%kg_k_gather,(sz1,sz2)) 1664 end if 1665 output%kg_k_gather(:,:)=reshape(buffer_int(ipck+1:ipck+nsize),(/sz1,sz2/)) 1666 ipck=ipck+nsize 1667 end if 1668 #else 1669 if (allocated(output%kg_k_gather)) then 1670 ABI_FREE(output%kg_k_gather) 1671 end if 1672 if (size1_kg_k_gather*size2_kg_k_gather>0) then 1673 nsize=size1_kg_k_gather*size2_kg_k_gather 1674 sz1=size1_kg_k_gather;sz2=size2_kg_k_gather 1675 ABI_MALLOC(output%kg_k_gather,(sz1,sz2)) 1676 output%kg_k_gather(:,:)=reshape(buffer_int(ipck+1:ipck+nsize),(/sz1,sz2/)) 1677 ipck=ipck+nsize 1678 end if 1679 #endif 1680 1681 if (allocated(output%recvcounts)) then 1682 ABI_FREE(output%recvcounts) 1683 end if 1684 if (size_recvcounts>0) then 1685 ABI_MALLOC(output%recvcounts,(size_recvcounts)) 1686 output%recvcounts(:)=buffer_int(ipck+1:ipck+size_recvcounts) 1687 ipck=ipck+size_recvcounts 1688 end if 1689 if (allocated(output%sendcounts)) then 1690 ABI_FREE(output%sendcounts) 1691 end if 1692 if (size_sendcounts>0) then 1693 ABI_MALLOC(output%sendcounts,(size_sendcounts)) 1694 output%sendcounts(:)=buffer_int(ipck+1:ipck+size_sendcounts) 1695 ipck=ipck+size_sendcounts 1696 end if 1697 if (allocated(output%rdispls)) then 1698 ABI_FREE(output%rdispls) 1699 end if 1700 if (size_rdispls>0) then 1701 ABI_MALLOC(output%rdispls,(size_rdispls)) 1702 output%rdispls(:)=buffer_int(ipck+1:ipck+size_rdispls) 1703 ipck=ipck+size_rdispls 1704 end if 1705 if (allocated(output%sdispls)) then 1706 ABI_FREE(output%sdispls) 1707 end if 1708 if (size_sdispls>0) then 1709 ABI_MALLOC(output%sdispls,(size_sdispls)) 1710 output%sdispls(:)=buffer_int(ipck+1:ipck+size_sdispls) 1711 ipck=ipck+size_sdispls 1712 end if 1713 if (allocated(output%gbound)) then 1714 ABI_FREE(output%gbound) 1715 end if 1716 if (size1_gbound*size2_gbound>0) then 1717 nsize=size1_gbound*size2_gbound 1718 sz1=size1_gbound;sz2=size2_gbound 1719 ABI_MALLOC(output%gbound,(sz1,sz2)) 1720 output%gbound(:,:)=reshape(buffer_int(ipck+1:ipck+nsize),(/sz1,sz2/)) 1721 ipck=ipck+nsize 1722 end if 1723 if (allocated(output%kg_k_gather_sym)) then 1724 ABI_FREE(output%kg_k_gather_sym) 1725 end if 1726 if (size1_kg_k_gather_sym*size2_kg_k_gather_sym>0) then 1727 nsize=size1_kg_k_gather_sym*size2_kg_k_gather_sym 1728 sz1=size1_kg_k_gather_sym;sz2=size2_kg_k_gather_sym 1729 ABI_MALLOC(output%kg_k_gather_sym,(sz1,sz2)) 1730 output%kg_k_gather_sym(:,:)=reshape(buffer_int(ipck+1:ipck+nsize),(/sz1,sz2/)) 1731 ipck=ipck+nsize 1732 end if 1733 if (allocated(output%rdispls_sym)) then 1734 ABI_FREE(output%rdispls_sym) 1735 end if 1736 if (size_rdispls_sym>0) then 1737 ABI_MALLOC(output%rdispls_sym,(size_rdispls_sym)) 1738 output%rdispls_sym(:)=buffer_int(ipck+1:ipck+size_rdispls_sym) 1739 ipck=ipck+size_rdispls_sym 1740 end if 1741 if (allocated(output%sdispls_sym)) then 1742 ABI_FREE(output%sdispls_sym) 1743 end if 1744 if (size_sdispls_sym>0) then 1745 ABI_MALLOC(output%sdispls_sym,(size_sdispls_sym)) 1746 output%sdispls_sym(:)=buffer_int(ipck+1:ipck+size_sdispls_sym) 1747 ipck=ipck+size_sdispls_sym 1748 end if 1749 if (allocated(output%recvcounts_sym)) then 1750 ABI_FREE(output%recvcounts_sym) 1751 end if 1752 if (size_recvcounts_sym>0) then 1753 ABI_MALLOC(output%recvcounts_sym,(size_recvcounts_sym)) 1754 output%recvcounts_sym(:)=buffer_int(ipck+1:ipck+size_recvcounts_sym) 1755 ipck=ipck+size_recvcounts_sym 1756 end if 1757 if (allocated(output%recvcounts_sym_tot)) then 1758 ABI_FREE(output%recvcounts_sym_tot) 1759 end if 1760 if (size_recvcounts_sym_tot>0) then 1761 ABI_MALLOC(output%recvcounts_sym_tot,(size_recvcounts_sym_tot)) 1762 output%recvcounts_sym_tot(:)=buffer_int(ipck+1:ipck+size_recvcounts_sym_tot) 1763 ipck=ipck+size_recvcounts_sym_tot 1764 end if 1765 if (allocated(output%sendcounts_sym)) then 1766 ABI_FREE(output%sendcounts_sym) 1767 end if 1768 if (size_sendcounts_sym>0) then 1769 ABI_MALLOC(output%sendcounts_sym,(size_sendcounts_sym)) 1770 output%sendcounts_sym(:)=buffer_int(ipck+1:ipck+size_sendcounts_sym) 1771 ipck=ipck+size_sendcounts_sym 1772 end if 1773 if (allocated(output%sendcounts_sym_all)) then 1774 ABI_FREE(output%sendcounts_sym_all) 1775 end if 1776 if (size_sendcounts_sym_all>0) then 1777 ABI_MALLOC(output%sendcounts_sym_all,(size_sendcounts_sym_all)) 1778 output%sendcounts_sym_all(:)=buffer_int(ipck+1:ipck+size_sendcounts_sym_all) 1779 ipck=ipck+size_sendcounts_sym_all 1780 end if 1781 if (allocated(output%tab_proc)) then 1782 ABI_FREE(output%tab_proc) 1783 end if 1784 if (size_tab_proc>0) then 1785 ABI_MALLOC(output%tab_proc,(size_tab_proc)) 1786 output%tab_proc(:)=buffer_int(ipck+1:ipck+size_tab_proc) 1787 ipck=ipck+size_tab_proc 1788 end if 1789 if (allocated(output%indices_pw_fft)) then 1790 ABI_FREE(output%indices_pw_fft) 1791 end if 1792 if (size_indices_pw_fft>0) then 1793 ABI_MALLOC(output%indices_pw_fft,(size_indices_pw_fft)) 1794 output%indices_pw_fft(:)=buffer_int(ipck+1:ipck+size_indices_pw_fft) 1795 ipck=ipck+size_indices_pw_fft 1796 end if 1797 if (allocated(output%sendcount_fft)) then 1798 ABI_FREE(output%sendcount_fft) 1799 end if 1800 if (size_sendcount_fft>0) then 1801 ABI_MALLOC(output%sendcount_fft,(size_sendcount_fft)) 1802 output%sendcount_fft(:)=buffer_int(ipck+1:ipck+size_sendcount_fft) 1803 ipck=ipck+size_sendcount_fft 1804 end if 1805 if (allocated(output%senddisp_fft)) then 1806 ABI_FREE(output%senddisp_fft) 1807 end if 1808 if (size_senddisp_fft>0) then 1809 ABI_MALLOC(output%senddisp_fft,(size_senddisp_fft)) 1810 output%senddisp_fft(:)=buffer_int(ipck+1:ipck+size_senddisp_fft) 1811 ipck=ipck+size_senddisp_fft 1812 end if 1813 if (allocated(output%recvcount_fft)) then 1814 ABI_FREE(output%recvcount_fft) 1815 end if 1816 if (size_recvcount_fft>0) then 1817 ABI_MALLOC(output%recvcount_fft,(size_recvcount_fft)) 1818 output%recvcount_fft(:)=buffer_int(ipck+1:ipck+size_recvcount_fft) 1819 ipck=ipck+size_recvcount_fft 1820 end if 1821 if (allocated(output%recvdisp_fft)) then 1822 ABI_FREE(output%recvdisp_fft) 1823 end if 1824 if (size_recvdisp_fft>0) then 1825 ABI_MALLOC(output%recvdisp_fft,(size_recvdisp_fft)) 1826 output%recvdisp_fft(:)=buffer_int(ipck+1:ipck+size_recvdisp_fft) 1827 ipck=ipck+size_recvdisp_fft 1828 end if 1829 if (allocated(output%kg_k_fft)) then 1830 ABI_FREE(output%kg_k_fft) 1831 end if 1832 if (size1_kg_k_fft*size2_kg_k_fft>0) then 1833 nsize=size1_kg_k_fft*size2_kg_k_fft 1834 sz1=size1_kg_k_fft;sz2=size2_kg_k_fft 1835 ABI_MALLOC(output%kg_k_fft,(sz1,sz2)) 1836 output%kg_k_fft(:,:)=reshape(buffer_int(ipck+1:ipck+nsize),(/sz1,sz2/)) 1837 ipck=ipck+nsize 1838 end if 1839 ABI_FREE(buffer_int) 1840 end if 1841 1842 !=== Receive and unpack real data ==== 1843 if (size_dp>0) then 1844 ipck=0 1845 ABI_MALLOC(buffer_dp,(size_dp)) 1846 call xmpi_recv(buffer_dp,sender,3*tag,spaceComm,ierr) 1847 if (allocated(output%ffnl_gather)) then 1848 ABI_FREE(output%ffnl_gather) 1849 end if 1850 if (size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather>0) then 1851 nsize=size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather 1852 sz1=size1_ffnl_gather;sz2=size2_ffnl_gather;sz3=size3_ffnl_gather;sz4=size4_ffnl_gather 1853 ABI_MALLOC(output%ffnl_gather,(sz1,sz2,sz3,sz4)) 1854 output%ffnl_gather(:,:,:,:)=reshape(buffer_dp(ipck+1:ipck+nsize),(/sz1,sz2,sz3,sz4/)) 1855 ipck=ipck+nsize 1856 end if 1857 1858 #if defined HAVE_GPU && defined HAVE_YAKL 1859 if (associated(output%kinpw_gather)) then 1860 if(output%gpu_option==ABI_GPU_KOKKOS) then 1861 ABI_FREE_MANAGED(output%kinpw_gather) 1862 else 1863 ABI_FREE(output%kinpw_gather) 1864 endif 1865 end if 1866 if (size_kinpw_gather>0) then 1867 if(output%gpu_option==ABI_GPU_KOKKOS) then 1868 ABI_MALLOC_MANAGED(output%kinpw_gather,(/size_kinpw_gather/)) 1869 else 1870 ABI_MALLOC(output%kinpw_gather,(size_kinpw_gather)) 1871 endif 1872 output%kinpw_gather(:)=buffer_dp(ipck+1:ipck+size_kinpw_gather) 1873 ipck=ipck+size_kinpw_gather 1874 end if 1875 #else 1876 if (allocated(output%kinpw_gather)) then 1877 ABI_FREE(output%kinpw_gather) 1878 end if 1879 if (size_kinpw_gather>0) then 1880 ABI_MALLOC(output%kinpw_gather,(size_kinpw_gather)) 1881 output%kinpw_gather(:)=buffer_dp(ipck+1:ipck+size_kinpw_gather) 1882 ipck=ipck+size_kinpw_gather 1883 end if 1884 #endif 1885 1886 if (allocated(output%ph3d_gather)) then 1887 ABI_FREE(output%ph3d_gather) 1888 end if 1889 if (size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather>0) then 1890 nsize=size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather 1891 sz1=size1_ph3d_gather;sz2=size2_ph3d_gather;sz3=size3_ph3d_gather 1892 ABI_MALLOC(output%ph3d_gather,(sz1,sz2,sz3)) 1893 output%ph3d_gather(:,:,:)=reshape(buffer_dp(ipck+1:ipck+nsize),(/sz1,sz2,sz3/)) 1894 ipck=ipck+nsize 1895 end if 1896 if (allocated(output%kpg_k_gather)) then 1897 ABI_FREE(output%kpg_k_gather) 1898 end if 1899 if (size1_kpg_k_gather*size2_kpg_k_gather>0) then 1900 nsize=size1_kpg_k_gather*size2_kpg_k_gather 1901 sz1=size1_kpg_k_gather;sz2=size2_kpg_k_gather 1902 ABI_MALLOC(output%kpg_k_gather,(sz1,sz2)) 1903 output%kpg_k_gather(:,:)=reshape(buffer_dp(ipck+1:ipck+nsize),(/sz1,sz2/)) 1904 ipck=ipck+nsize 1905 end if 1906 ABI_FREE(buffer_dp) 1907 end if 1908 1909 end subroutine bandfft_kpt_mpi_recv
ABINIT/bandfft_kpt_mpi_send [ Functions ]
NAME
bandfft_kpt_mpi_send
FUNCTION
Send a bandfft_kpt_type inside a MPI communicator
INPUTS
input=The datatype to be transmitted receiver=ID of the receiver in spaceComm tag=message tag spaceComm=MPI Communicator [profile]=(character, optional, default=full) Define the part of the datastructure to be sent; possible values: 'full'= the entire datastructure 'fourwf'= only arrays needed by the fourwf routine
OUTPUT
ierr=Error status
SOURCE
1244 subroutine bandfft_kpt_mpi_send(input,receiver,tag,spaceComm,ierr,profile) 1245 1246 !Arguments ------------------------------------ 1247 !scalars 1248 integer,intent(in) :: receiver,spaceComm,tag 1249 integer,intent(out) :: ierr 1250 character(len=*),optional,intent(in) :: profile 1251 !arrays 1252 type(bandfft_kpt_type),intent(in) :: input 1253 1254 !Local variables------------------------------- 1255 !scalars 1256 integer :: ipck,nsize,size_dp,size_int 1257 integer :: size1_kg_k_gather=0,size2_kg_k_gather=0,size_recvcounts=0 1258 integer :: size_sendcounts=0,size_rdispls=0,size_sdispls=0,size1_gbound=0,size2_gbound 1259 integer :: size1_kg_k_gather_sym=0,size2_kg_k_gather_sym=0,size_rdispls_sym=0 1260 integer :: size_sdispls_sym=0,size_recvcounts_sym=0,size_recvcounts_sym_tot=0 1261 integer :: size_sendcounts_sym=0,size_sendcounts_sym_all=0,size_tab_proc=0 1262 integer :: size_indices_pw_fft=0,size_sendcount_fft=0,size_senddisp_fft=0 1263 integer :: size_recvcount_fft=0,size_recvdisp_fft=0,size1_kg_k_fft=0,size2_kg_k_fft=0 1264 integer :: size1_ffnl_gather=0,size2_ffnl_gather=0,size3_ffnl_gather=0,size4_ffnl_gather=0 1265 integer :: size_kinpw_gather=0,size1_ph3d_gather=0,size2_ph3d_gather=0,size3_ph3d_gather=0 1266 integer :: size1_kpg_k_gather=0,size2_kpg_k_gather=0 1267 logical :: fourwf,full 1268 !arrays 1269 integer,allocatable :: buffer_int(:) 1270 real(dp),allocatable :: buffer_dp(:) 1271 1272 ! ************************************************************************* 1273 1274 if (xmpi_comm_size(spaceComm)<=1) return 1275 1276 fourwf=.false.;if (present(profile)) fourwf=(trim(profile)=='fourwf') 1277 full=(.not.fourwf) 1278 1279 !=== Store sizes ==== 1280 if (fourwf.or.full) then 1281 1282 #if defined HAVE_GPU && defined HAVE_YAKL 1283 if (associated(input%kg_k_gather)) size1_kg_k_gather=size(input%kg_k_gather,1) 1284 if (associated(input%kg_k_gather)) size2_kg_k_gather=size(input%kg_k_gather,2) 1285 #else 1286 if (allocated(input%kg_k_gather)) size1_kg_k_gather=size(input%kg_k_gather,1) 1287 if (allocated(input%kg_k_gather)) size2_kg_k_gather=size(input%kg_k_gather,2) 1288 #endif 1289 if (input%flag1_is_allocated==1) then 1290 if (allocated(input%recvcounts)) size_recvcounts=size(input%recvcounts) 1291 if (allocated(input%sendcounts)) size_sendcounts=size(input%sendcounts) 1292 if (allocated(input%rdispls)) size_rdispls=size(input%rdispls) 1293 if (allocated(input%sdispls)) size_sdispls=size(input%sdispls) 1294 if (allocated(input%gbound)) size1_gbound=size(input%gbound,1) 1295 if (allocated(input%gbound)) size2_gbound=size(input%gbound,2) 1296 end if 1297 if (input%flag3_is_allocated==1.and.input%istwf_k>1) then 1298 if (allocated(input%kg_k_gather_sym)) size1_kg_k_gather_sym=size(input%kg_k_gather_sym,1) 1299 if (allocated(input%kg_k_gather_sym)) size2_kg_k_gather_sym=size(input%kg_k_gather_sym,2) 1300 if (allocated(input%rdispls_sym)) size_rdispls_sym=size(input%rdispls_sym) 1301 if (allocated(input%sdispls_sym)) size_sdispls_sym=size(input%sdispls_sym) 1302 if (allocated(input%recvcounts_sym)) size_recvcounts_sym=size(input%recvcounts_sym) 1303 if (allocated(input%recvcounts_sym_tot)) size_recvcounts_sym_tot=size(input%recvcounts_sym_tot) 1304 if (allocated(input%sendcounts_sym)) size_sendcounts_sym=size(input%sendcounts_sym) 1305 if (allocated(input%sendcounts_sym_all)) size_sendcounts_sym_all=size(input%sendcounts_sym_all) 1306 if (allocated(input%tab_proc)) size_tab_proc=size(input%tab_proc) 1307 end if 1308 if (input%have_to_reequilibrate) then 1309 if (allocated(input%indices_pw_fft)) size_indices_pw_fft=size(input%indices_pw_fft) 1310 if (allocated(input%sendcount_fft)) size_sendcount_fft=size(input%sendcount_fft) 1311 if (allocated(input%senddisp_fft)) size_senddisp_fft=size(input%senddisp_fft) 1312 if (allocated(input%recvcount_fft)) size_recvcount_fft=size(input%recvcount_fft) 1313 if (allocated(input%recvdisp_fft)) size_recvdisp_fft=size(input%recvdisp_fft) 1314 if (allocated(input%kg_k_fft)) size1_kg_k_fft=size(input%kg_k_fft,1) 1315 if (allocated(input%kg_k_fft)) size2_kg_k_fft=size(input%kg_k_fft,2) 1316 end if 1317 end if 1318 if (input%flag2_is_allocated==1.and.full) then 1319 if (allocated(input%ffnl_gather)) size1_ffnl_gather=size(input%ffnl_gather,1) 1320 if (allocated(input%ffnl_gather)) size2_ffnl_gather=size(input%ffnl_gather,2) 1321 if (allocated(input%ffnl_gather)) size3_ffnl_gather=size(input%ffnl_gather,3) 1322 if (allocated(input%ffnl_gather)) size4_ffnl_gather=size(input%ffnl_gather,4) 1323 1324 #if defined HAVE_GPU && defined HAVE_YAKL 1325 if (associated(input%kinpw_gather)) size_kinpw_gather=size(input%kinpw_gather) 1326 #else 1327 if (allocated(input%kinpw_gather)) size_kinpw_gather=size(input%kinpw_gather) 1328 #endif 1329 1330 if (allocated(input%ph3d_gather)) size1_ph3d_gather=size(input%ph3d_gather,1) 1331 if (allocated(input%ph3d_gather)) size2_ph3d_gather=size(input%ph3d_gather,2) 1332 if (allocated(input%ph3d_gather)) size3_ph3d_gather=size(input%ph3d_gather,3) 1333 if (allocated(input%kpg_k_gather)) size1_kpg_k_gather=size(input%kpg_k_gather,1) 1334 if (allocated(input%kpg_k_gather)) size2_kpg_k_gather=size(input%kpg_k_gather,2) 1335 end if 1336 1337 !=== Compute amount of data to transmit ==== 1338 size_int=size1_kg_k_gather*size2_kg_k_gather & 1339 & +size_recvcounts+size_sendcounts+size_rdispls+size_sdispls & 1340 & +size1_gbound*size2_gbound & 1341 & +size1_kg_k_gather_sym*size2_kg_k_gather_sym & 1342 & +size_rdispls_sym+size_sdispls_sym & 1343 & +size_recvcounts_sym+size_recvcounts_sym_tot & 1344 & +size_sendcounts_sym+size_sendcounts_sym_all+size_tab_proc & 1345 & +size_indices_pw_fft+size_sendcount_fft+size_senddisp_fft & 1346 & +size_recvcount_fft+size_recvdisp_fft & 1347 & +size1_kg_k_fft*size2_kg_k_fft 1348 size_dp=size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather & 1349 & +size_kinpw_gather & 1350 & +size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather & 1351 & +size1_kpg_k_gather*size2_kpg_k_gather 1352 1353 !=== Send flags and array sizes ==== 1354 nsize=45 1355 ABI_MALLOC(buffer_int,(nsize)) 1356 buffer_int(1)=input%flag1_is_allocated 1357 buffer_int(2)=input%flag2_is_allocated 1358 buffer_int(3)=input%flag3_is_allocated 1359 buffer_int(4)=0;if (input%have_to_reequilibrate) buffer_int(4)=1 1360 buffer_int(5)=input%istwf_k 1361 buffer_int(6)=input%npw_tot 1362 buffer_int(7)=input%npw_fft 1363 buffer_int(8)=input%ndatarecv 1364 buffer_int(9)=input%idatarecv0 1365 buffer_int(10)=input%ndatarecv_tot 1366 buffer_int(11)=input%ndatasend_sym 1367 buffer_int(12)=size_recvcounts 1368 buffer_int(13)=size_sendcounts 1369 buffer_int(14)=size_rdispls 1370 buffer_int(15)=size_sdispls 1371 buffer_int(16)=size1_gbound 1372 buffer_int(17)=size2_gbound 1373 buffer_int(18)=size1_ffnl_gather 1374 buffer_int(19)=size2_ffnl_gather 1375 buffer_int(20)=size3_ffnl_gather 1376 buffer_int(21)=size4_ffnl_gather 1377 buffer_int(22)=size_kinpw_gather 1378 buffer_int(23)=size1_ph3d_gather 1379 buffer_int(24)=size2_ph3d_gather 1380 buffer_int(25)=size3_ph3d_gather 1381 buffer_int(26)=size1_kpg_k_gather 1382 buffer_int(27)=size2_kpg_k_gather 1383 buffer_int(28)=size_indices_pw_fft 1384 buffer_int(29)=size_sendcount_fft 1385 buffer_int(30)=size_senddisp_fft 1386 buffer_int(31)=size_recvcount_fft 1387 buffer_int(32)=size_recvdisp_fft 1388 buffer_int(33)=size1_kg_k_fft 1389 buffer_int(34)=size2_kg_k_fft 1390 buffer_int(35)=size1_kg_k_gather 1391 buffer_int(36)=size2_kg_k_gather 1392 buffer_int(37)=size1_kg_k_gather_sym 1393 buffer_int(38)=size2_kg_k_gather_sym 1394 buffer_int(39)=size_rdispls_sym 1395 buffer_int(40)=size_sdispls_sym 1396 buffer_int(41)=size_recvcounts_sym 1397 buffer_int(42)=size_recvcounts_sym_tot 1398 buffer_int(43)=size_sendcounts_sym 1399 buffer_int(44)=size_sendcounts_sym_all 1400 buffer_int(45)=size_tab_proc 1401 call xmpi_send(buffer_int,receiver,3*tag-2,spaceComm,ierr) 1402 ABI_FREE(buffer_int) 1403 1404 !=== Pack and send integer data ==== 1405 if (size_int>0) then 1406 ipck=0 1407 ABI_MALLOC(buffer_int,(size_int)) 1408 if (size1_kg_k_gather*size2_kg_k_gather>0) then 1409 nsize=size1_kg_k_gather*size2_kg_k_gather 1410 buffer_int(ipck+1:ipck+nsize)=reshape(input%kg_k_gather(:,:),(/nsize/)) 1411 ipck=ipck+nsize 1412 end if 1413 if (size_recvcounts>0) then 1414 buffer_int(ipck+1:ipck+size_recvcounts)=input%recvcounts(:) 1415 ipck=ipck+size_recvcounts 1416 end if 1417 if (size_sendcounts>0) then 1418 buffer_int(ipck+1:ipck+size_sendcounts)=input%sendcounts(:) 1419 ipck=ipck+size_sendcounts 1420 end if 1421 if (size_rdispls>0) then 1422 buffer_int(ipck+1:ipck+size_rdispls)=input%rdispls(:) 1423 ipck=ipck+size_rdispls 1424 end if 1425 if (size_sdispls>0) then 1426 buffer_int(ipck+1:ipck+size_sdispls)=input%sdispls(:) 1427 ipck=ipck+size_sdispls 1428 end if 1429 if (size1_gbound*size2_gbound>0) then 1430 nsize=size1_gbound*size2_gbound 1431 buffer_int(ipck+1:ipck+nsize)=reshape(input%gbound(:,:),(/nsize/)) 1432 ipck=ipck+nsize 1433 end if 1434 if (size1_kg_k_gather_sym*size2_kg_k_gather_sym>0) then 1435 nsize=size1_kg_k_gather_sym*size2_kg_k_gather_sym 1436 buffer_int(ipck+1:ipck+nsize)=reshape(input%kg_k_gather_sym(:,:),(/nsize/)) 1437 ipck=ipck+nsize 1438 end if 1439 if (size_rdispls_sym>0) then 1440 buffer_int(ipck+1:ipck+size_rdispls_sym)=input%rdispls_sym(:) 1441 ipck=ipck+size_rdispls_sym 1442 end if 1443 if (size_sdispls_sym>0) then 1444 buffer_int(ipck+1:ipck+size_sdispls_sym)=input%sdispls_sym(:) 1445 ipck=ipck+size_sdispls_sym 1446 end if 1447 if (size_recvcounts_sym>0) then 1448 buffer_int(ipck+1:ipck+size_recvcounts_sym)=input%recvcounts_sym(:) 1449 ipck=ipck+size_recvcounts_sym 1450 end if 1451 if (size_recvcounts_sym_tot>0) then 1452 buffer_int(ipck+1:ipck+size_recvcounts_sym_tot)=input%recvcounts_sym_tot(:) 1453 ipck=ipck+size_recvcounts_sym_tot 1454 end if 1455 if (size_sendcounts_sym>0) then 1456 buffer_int(ipck+1:ipck+size_sendcounts_sym)=input%sendcounts_sym(:) 1457 ipck=ipck+size_sendcounts_sym 1458 end if 1459 if (size_sendcounts_sym_all>0) then 1460 buffer_int(ipck+1:ipck+size_sendcounts_sym_all)=input%sendcounts_sym_all(:) 1461 ipck=ipck+size_sendcounts_sym_all 1462 end if 1463 if (size_tab_proc>0) then 1464 buffer_int(ipck+1:ipck+size_tab_proc)=input%tab_proc(:) 1465 ipck=ipck+size_tab_proc 1466 end if 1467 if (size_indices_pw_fft>0) then 1468 buffer_int(ipck+1:ipck+size_indices_pw_fft)=input%indices_pw_fft(:) 1469 ipck=ipck+size_indices_pw_fft 1470 end if 1471 if (size_sendcount_fft>0) then 1472 buffer_int(ipck+1:ipck+size_sendcount_fft)=input%sendcount_fft(:) 1473 ipck=ipck+size_sendcount_fft 1474 end if 1475 if (size_senddisp_fft>0) then 1476 buffer_int(ipck+1:ipck+size_senddisp_fft)=input%senddisp_fft(:) 1477 ipck=ipck+size_senddisp_fft 1478 end if 1479 if (size_recvcount_fft>0) then 1480 buffer_int(ipck+1:ipck+size_recvcount_fft)=input%recvcount_fft(:) 1481 ipck=ipck+size_recvcount_fft 1482 end if 1483 if (size_recvdisp_fft>0) then 1484 buffer_int(ipck+1:ipck+size_recvdisp_fft)=input%recvdisp_fft(:) 1485 ipck=ipck+size_recvdisp_fft 1486 end if 1487 if (size1_kg_k_fft*size2_kg_k_fft>0) then 1488 nsize=size1_kg_k_fft*size2_kg_k_fft 1489 buffer_int(ipck+1:ipck+nsize)=reshape(input%kg_k_fft(:,:),(/nsize/)) 1490 ipck=ipck+nsize 1491 end if 1492 call xmpi_send(buffer_int,receiver,3*tag-1,spaceComm,ierr) 1493 ABI_FREE(buffer_int) 1494 end if 1495 1496 !=== Pack and send real data ==== 1497 if (size_dp>0) then 1498 ipck=0 1499 ABI_MALLOC(buffer_dp,(size_dp)) 1500 if (size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather>0) then 1501 nsize=size1_ffnl_gather*size2_ffnl_gather*size3_ffnl_gather*size4_ffnl_gather 1502 buffer_dp(ipck+1:ipck+nsize)=reshape(input%ffnl_gather(:,:,:,:),(/nsize/)) 1503 ipck=ipck+nsize 1504 end if 1505 if (size_kinpw_gather>0) then 1506 buffer_dp(ipck+1:ipck+size_kinpw_gather)=input%kinpw_gather(:) 1507 ipck=ipck+size_kinpw_gather 1508 end if 1509 if (size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather>0) then 1510 nsize=size1_ph3d_gather*size2_ph3d_gather*size3_ph3d_gather 1511 buffer_dp(ipck+1:ipck+nsize)=reshape(input%ph3d_gather(:,:,:),(/nsize/)) 1512 ipck=ipck+nsize 1513 end if 1514 if (size1_kpg_k_gather*size2_kpg_k_gather>0) then 1515 nsize=size1_kpg_k_gather*size2_kpg_k_gather 1516 buffer_dp(ipck+1:ipck+nsize)=reshape(input%kpg_k_gather(:,:),(/nsize/)) 1517 ipck=ipck+nsize 1518 end if 1519 call xmpi_send(buffer_dp,receiver,3*tag,spaceComm,ierr) 1520 ABI_FREE(buffer_dp) 1521 end if 1522 1523 end subroutine bandfft_kpt_mpi_send
ABINIT/m_bandfft_kpt [ Modules ]
NAME
m_bandfft_kpt
FUNCTION
This module provides the definition of the bandfft_kpt_type used for kgb parallelization.
COPYRIGHT
Copyright (C) 2011-2024 ABINIT group (FJ, FB, 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
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 MODULE m_bandfft_kpt 24 25 use, intrinsic :: iso_c_binding, only : c_int32_t, c_double 26 27 use defs_basis 28 use m_abicore 29 use m_errors 30 use m_xmpi 31 32 use defs_abitypes, only : MPI_type 33 use m_time, only : timab 34 use m_kg, only : mkkpg 35 use m_fftcore, only : sphereboundary 36 use m_mpinfo, only : proc_distrb_cycle 37 use m_hamiltonian, only : gs_hamiltonian_type 38 39 #if defined HAVE_YAKL 40 use gator_mod 41 #endif 42 43 implicit none 44 45 private 46 public :: bandfft_kpt_init1 47 public :: bandfft_kpt_init2 48 public :: bandfft_kpt_reset 49 public :: bandfft_kpt_destroy 50 public :: bandfft_kpt_destroy_array 51 public :: bandfft_kpt_copy 52 public :: bandfft_kpt_mpi_send 53 public :: bandfft_kpt_mpi_recv 54 public :: bandfft_kpt_savetabs 55 public :: bandfft_kpt_restoretabs 56 public :: bandfft_kpt_set_ikpt 57 public :: bandfft_kpt_get_ikpt 58 public :: prep_bandfft_tabs
ABINIT/prep_bandfft_tabs [ Functions ]
NAME
prep_bandfft_tabs
FUNCTION
This routine transpose various tabs needed in bandfft parallelization
INPUTS
ikpt=index of the k-point gs_hamk <type(gs_hamiltonian_type)>=all data for the Hamiltonian at k mkmem=number of k points which can fit in memory; set to 0 if use disk
SIDE EFFECTS
bandfft_kpt tabs (defined in m_bandfft_kpt module)
SOURCE
2199 subroutine prep_bandfft_tabs(gs_hamk,ikpt,mkmem,mpi_enreg) 2200 2201 !Arguments ------------------------------- 2202 integer,intent(in) :: ikpt,mkmem 2203 type(gs_hamiltonian_type),intent(inout) :: gs_hamk 2204 type(MPI_type),intent(inout) :: mpi_enreg 2205 2206 !Local variables------------------------------- 2207 integer :: dimffnl,ierr,ikpt_this_proc,ipw,lmnmax,matblk,ndatarecv,nkpg,npw_k,ntypat,spaceComm 2208 logical :: tabs_allocated 2209 real(dp) :: tsec(2) 2210 character(len=500) :: message 2211 integer, allocatable :: recvcounts(:),rdispls(:) 2212 integer, allocatable :: recvcountsloc(:),rdisplsloc(:) 2213 real(dp),allocatable :: ffnl_gather(:,:,:,:),ffnl_little(:,:,:,:),ffnl_little_gather(:,:,:,:) 2214 real(dp),allocatable :: kinpw_gather(:),kpg_k_gather(:,:) 2215 real(dp),allocatable :: ph3d_gather(:,:,:),ph3d_little(:,:,:),ph3d_little_gather(:,:,:) 2216 2217 ! ********************************************************************* 2218 2219 call timab(575,1,tsec) 2220 2221 ikpt_this_proc=mpi_enreg%my_kpttab(ikpt) 2222 tabs_allocated=((bandfft_kpt(ikpt_this_proc)%flag1_is_allocated==1).and.& 2223 & (ikpt_this_proc <= mkmem).and.(ikpt_this_proc/=0)) 2224 2225 if (.not.tabs_allocated) then 2226 message = ' the bandfft tabs are not allocated !' 2227 ABI_BUG(message) 2228 end if 2229 2230 ntypat=gs_hamk%ntypat 2231 lmnmax=gs_hamk%lmnmax 2232 matblk=gs_hamk%matblk 2233 npw_k=gs_hamk%npw_k 2234 dimffnl=size(gs_hamk%ffnl_k,2) 2235 nkpg=size(gs_hamk%kpg_k,2) 2236 2237 ABI_MALLOC(rdispls,(mpi_enreg%nproc_band)) 2238 ABI_MALLOC(recvcounts,(mpi_enreg%nproc_band)) 2239 ABI_MALLOC(rdisplsloc,(mpi_enreg%nproc_band)) 2240 ABI_MALLOC(recvcountsloc,(mpi_enreg%nproc_band)) 2241 2242 spaceComm =mpi_enreg%comm_band 2243 ndatarecv =bandfft_kpt(ikpt_this_proc)%ndatarecv 2244 rdispls(:) =bandfft_kpt(ikpt_this_proc)%rdispls(:) 2245 recvcounts(:)=bandfft_kpt(ikpt_this_proc)%recvcounts(:) 2246 2247 !---- Process FFNL 2248 if (associated(gs_hamk%ffnl_k)) then 2249 ABI_MALLOC(ffnl_gather,(ndatarecv,dimffnl,lmnmax,ntypat)) 2250 ABI_MALLOC(ffnl_little,(dimffnl,lmnmax,ntypat,npw_k)) 2251 ABI_MALLOC(ffnl_little_gather,(dimffnl,lmnmax,ntypat,ndatarecv)) 2252 do ipw=1,npw_k 2253 ffnl_little(:,:,:,ipw)=gs_hamk%ffnl_k(ipw,:,:,:) 2254 end do 2255 recvcountsloc(:)=recvcounts(:)*dimffnl*lmnmax*ntypat 2256 rdisplsloc(:)=rdispls(:)*dimffnl*lmnmax*ntypat 2257 call xmpi_allgatherv(ffnl_little,npw_k*dimffnl*lmnmax*ntypat,ffnl_little_gather,& 2258 & recvcountsloc(:),rdisplsloc,spaceComm,ierr) 2259 do ipw=1,ndatarecv 2260 ffnl_gather(ipw,:,:,:)=ffnl_little_gather(:,:,:,ipw) 2261 end do 2262 ABI_FREE(ffnl_little) 2263 ABI_FREE(ffnl_little_gather) 2264 else 2265 ABI_MALLOC(ffnl_gather,(0,0,0,0)) 2266 end if 2267 2268 !---- Process PH3D 2269 if (associated(gs_hamk%ph3d_k)) then 2270 ABI_MALLOC(ph3d_gather,(2,ndatarecv,matblk)) 2271 ABI_MALLOC(ph3d_little,(2,matblk,npw_k)) 2272 ABI_MALLOC(ph3d_little_gather,(2,matblk,ndatarecv)) 2273 recvcountsloc(:)=recvcounts(:)*2*matblk 2274 rdisplsloc(:)=rdispls(:)*2*matblk 2275 do ipw=1,npw_k 2276 ph3d_little(:,:,ipw)=gs_hamk%ph3d_k(:,ipw,:) 2277 end do 2278 call xmpi_allgatherv(ph3d_little,npw_k*2*matblk,ph3d_little_gather,& 2279 & recvcountsloc(:),rdisplsloc,spaceComm,ierr) 2280 do ipw=1,ndatarecv 2281 ph3d_gather(:,ipw,:)=ph3d_little_gather(:,:,ipw) 2282 end do 2283 ABI_FREE(ph3d_little_gather) 2284 ABI_FREE(ph3d_little) 2285 else 2286 ABI_MALLOC(ph3d_gather,(0,0,0)) 2287 end if 2288 2289 !---- Process KPG_K 2290 if (associated(gs_hamk%kpg_k)) then 2291 ABI_MALLOC(kpg_k_gather,(ndatarecv,nkpg)) 2292 if (nkpg>0) then 2293 call mkkpg(bandfft_kpt(ikpt_this_proc)%kg_k_gather,kpg_k_gather,gs_hamk%kpt_k,nkpg,ndatarecv) 2294 end if 2295 else 2296 ABI_MALLOC(kpg_k_gather,(0,0)) 2297 end if 2298 2299 !---- Process KINPW 2300 if (associated(gs_hamk%kinpw_k)) then 2301 ABI_MALLOC(kinpw_gather,(ndatarecv)) 2302 recvcountsloc(:)=recvcounts(:) 2303 rdisplsloc(:)=rdispls(:) 2304 call xmpi_allgatherv(gs_hamk%kinpw_k,npw_k,kinpw_gather,recvcountsloc(:),rdisplsloc,spaceComm,ierr) 2305 else 2306 ABI_MALLOC(kinpw_gather,(0)) 2307 end if 2308 2309 ABI_FREE(recvcounts) 2310 ABI_FREE(rdispls) 2311 ABI_FREE(recvcountsloc) 2312 ABI_FREE(rdisplsloc) 2313 2314 call bandfft_kpt_init2(bandfft_kpt,dimffnl,ffnl_gather,ikpt_this_proc,kinpw_gather,kpg_k_gather,& 2315 & lmnmax,matblk,mkmem,ndatarecv,nkpg,ntypat,ph3d_gather) 2316 2317 ABI_FREE(ffnl_gather) 2318 ABI_FREE(ph3d_gather) 2319 ABI_FREE(kinpw_gather) 2320 ABI_FREE(kpg_k_gather) 2321 2322 !---- Store current kpt index 2323 call bandfft_kpt_set_ikpt(ikpt,mpi_enreg) 2324 2325 call timab(575,2,tsec) 2326 2327 end subroutine prep_bandfft_tabs
m_bandfft_kpt/bandfft_kpt_copy [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_copy
FUNCTION
Copy a bandfft_kpt datastructure into another
INPUTS
bandfft_kpt_in=<type(bandfft_kpt_type)>=input bandfft_kpt datastructure
OUTPUT
bandfft_kpt_out=<type(bandfft_kpt_type)>=output bandfft_kpt datastructure
SOURCE
1002 subroutine bandfft_kpt_copy(bandfft_kpt_in,bandfft_kpt_out,mpi_enreg1,opt_bandfft) 1003 1004 !Arguments ------------------------------------ 1005 !scalars 1006 integer,intent(in) :: opt_bandfft 1007 type(bandfft_kpt_type),pointer :: bandfft_kpt_in(:) 1008 type(bandfft_kpt_type),pointer :: bandfft_kpt_out(:) 1009 type(MPI_type),intent(inout) :: mpi_enreg1 1010 1011 !Local variables------------------------------- 1012 !scalars 1013 integer :: ikpt,isppol,jkpt,sz1,sz2,sz3,sz4 1014 1015 ! ********************************************************************* 1016 1017 !Optional pointers 1018 if (opt_bandfft==0) then 1019 nullify(bandfft_kpt_out) 1020 else if (opt_bandfft==1) then 1021 if (associated(bandfft_kpt_in)) then 1022 ABI_MALLOC(bandfft_kpt_out,(size(bandfft_kpt_in))) 1023 do isppol=1,size(mpi_enreg1%proc_distrb,3) 1024 do ikpt=1,size(mpi_enreg1%proc_distrb,1) 1025 sz1=size(mpi_enreg1%proc_distrb,2) 1026 if(proc_distrb_cycle(mpi_enreg1%proc_distrb,ikpt,1,sz1,isppol,mpi_enreg1%me_kpt)) then 1027 cycle 1028 end if 1029 jkpt=mpi_enreg1%my_kpttab(ikpt) 1030 ! if (allocated(bandfft_kpt_in(jkpt)%ind_kg_mpi_to_seq)) then 1031 ! sz1=size(bandfft_kpt_in(jkpt)%ind_kg_mpi_to_seq) 1032 ! ABI_MALLOC(bandfft_kpt_out(jkpt)%ind_kg_mpi_to_seq,(sz1)) 1033 ! bandfft_kpt_out(jkpt)%ind_kg_mpi_to_seq= & 1034 ! & bandfft_kpt_in(jkpt)%ind_kg_mpi_to_seq 1035 ! end if 1036 1037 bandfft_kpt_out(jkpt)%gpu_option=bandfft_kpt_in(jkpt)%gpu_option 1038 #if defined HAVE_GPU && defined HAVE_YAKL 1039 if (associated(bandfft_kpt_in(jkpt)%kg_k_gather)) then 1040 sz1=size(bandfft_kpt_in(jkpt)%kg_k_gather,1) 1041 sz2=size(bandfft_kpt_in(jkpt)%kg_k_gather,2) 1042 if(bandfft_kpt_in(jkpt)%gpu_option==ABI_GPU_KOKKOS) then 1043 ABI_MALLOC_MANAGED(bandfft_kpt_out(jkpt)%kg_k_gather,(/sz1,sz2/)) 1044 bandfft_kpt_out(jkpt)%kg_k_gather= & 1045 & bandfft_kpt_in(jkpt)%kg_k_gather 1046 else 1047 ABI_MALLOC(bandfft_kpt_out(jkpt)%kg_k_gather,(sz1,sz2)) 1048 bandfft_kpt_out(jkpt)%kg_k_gather= & 1049 & bandfft_kpt_in(jkpt)%kg_k_gather 1050 end if 1051 end if 1052 #else 1053 if (allocated(bandfft_kpt_in(jkpt)%kg_k_gather)) then 1054 sz1=size(bandfft_kpt_in(jkpt)%kg_k_gather,1) 1055 sz2=size(bandfft_kpt_in(jkpt)%kg_k_gather,2) 1056 ABI_MALLOC(bandfft_kpt_out(jkpt)%kg_k_gather,(sz1,sz2)) 1057 bandfft_kpt_out(jkpt)%kg_k_gather= & 1058 & bandfft_kpt_in(jkpt)%kg_k_gather 1059 end if 1060 #endif 1061 1062 bandfft_kpt_out(jkpt)%flag1_is_allocated=bandfft_kpt_in(jkpt)%flag1_is_allocated 1063 if (allocated(bandfft_kpt_in(jkpt)%gbound)) then 1064 sz1=size(bandfft_kpt_in(jkpt)%gbound,1) 1065 sz2=size(bandfft_kpt_in(jkpt)%gbound,2) 1066 ABI_MALLOC(bandfft_kpt_out(jkpt)%gbound,(sz1,sz2)) 1067 bandfft_kpt_out(jkpt)%gbound= & 1068 & bandfft_kpt_in(jkpt)%gbound 1069 end if 1070 if (allocated(bandfft_kpt_in(jkpt)%recvcounts)) then 1071 sz1=size(bandfft_kpt_in(jkpt)%recvcounts) 1072 ABI_MALLOC(bandfft_kpt_out(jkpt)%recvcounts,(sz1)) 1073 bandfft_kpt_out(jkpt)%recvcounts= & 1074 & bandfft_kpt_in(jkpt)%recvcounts 1075 end if 1076 if (allocated(bandfft_kpt_in(jkpt)%sendcounts)) then 1077 ABI_MALLOC(bandfft_kpt_out(jkpt)%sendcounts,(size(bandfft_kpt_in(jkpt)%sendcounts))) 1078 bandfft_kpt_out(jkpt)%sendcounts= & 1079 & bandfft_kpt_in(jkpt)%sendcounts 1080 end if 1081 if (allocated(bandfft_kpt_in(jkpt)%rdispls)) then 1082 ABI_MALLOC(bandfft_kpt_out(jkpt)%rdispls,(size(bandfft_kpt_in(jkpt)%rdispls))) 1083 bandfft_kpt_out(jkpt)%rdispls= & 1084 & bandfft_kpt_in(jkpt)%rdispls 1085 end if 1086 if (allocated(bandfft_kpt_in(jkpt)%sdispls)) then 1087 ABI_MALLOC(bandfft_kpt_out(jkpt)%sdispls,(size(bandfft_kpt_in(jkpt)%sdispls))) 1088 bandfft_kpt_out(jkpt)%sdispls= & 1089 & bandfft_kpt_in(jkpt)%sdispls 1090 end if 1091 bandfft_kpt_out(jkpt)%flag2_is_allocated=bandfft_kpt_in(jkpt)%flag2_is_allocated 1092 if (allocated(bandfft_kpt_in(jkpt)%ffnl_gather)) then 1093 sz1=size(bandfft_kpt_in(jkpt)%ffnl_gather,1) 1094 sz2=size(bandfft_kpt_in(jkpt)%ffnl_gather,2) 1095 sz3=size(bandfft_kpt_in(jkpt)%ffnl_gather,3) 1096 sz4=size(bandfft_kpt_in(jkpt)%ffnl_gather,4) 1097 ABI_MALLOC(bandfft_kpt_out(jkpt)%ffnl_gather,(sz1,sz2,sz3,sz4)) 1098 bandfft_kpt_out(jkpt)%ffnl_gather= & 1099 & bandfft_kpt_in(jkpt)%ffnl_gather 1100 end if 1101 1102 #if defined HAVE_GPU && defined HAVE_YAKL 1103 if (associated(bandfft_kpt_in(jkpt)%kinpw_gather)) then 1104 sz1=size(bandfft_kpt_in(jkpt)%kinpw_gather) 1105 if(bandfft_kpt_in(jkpt)%gpu_option==ABI_GPU_KOKKOS) then 1106 ABI_MALLOC_MANAGED(bandfft_kpt_out(jkpt)%kinpw_gather,(/sz1/)) 1107 bandfft_kpt_out(jkpt)%kinpw_gather= & 1108 & bandfft_kpt_in(jkpt)%kinpw_gather 1109 else 1110 ABI_MALLOC(bandfft_kpt_out(jkpt)%kinpw_gather,(sz1)) 1111 bandfft_kpt_out(jkpt)%kinpw_gather= & 1112 & bandfft_kpt_in(jkpt)%kinpw_gather 1113 end if 1114 end if 1115 #else 1116 if (allocated(bandfft_kpt_in(jkpt)%kinpw_gather)) then 1117 sz1=size(bandfft_kpt_in(jkpt)%kinpw_gather) 1118 ABI_MALLOC(bandfft_kpt_out(jkpt)%kinpw_gather,(sz1)) 1119 bandfft_kpt_out(jkpt)%kinpw_gather= & 1120 & bandfft_kpt_in(jkpt)%kinpw_gather 1121 end if 1122 #endif 1123 1124 if (allocated(bandfft_kpt_in(jkpt)%ph3d_gather)) then 1125 sz1=size(bandfft_kpt_in(jkpt)%ph3d_gather,1) 1126 sz2=size(bandfft_kpt_in(jkpt)%ph3d_gather,2) 1127 sz3=size(bandfft_kpt_in(jkpt)%ph3d_gather,3) 1128 ABI_MALLOC(bandfft_kpt_out(jkpt)%ph3d_gather,(sz1,sz2,sz3)) 1129 bandfft_kpt_out(jkpt)%ph3d_gather= & 1130 & bandfft_kpt_in(jkpt)%ph3d_gather 1131 end if 1132 if (allocated(bandfft_kpt_in(jkpt)%kpg_k_gather)) then 1133 sz1=size(bandfft_kpt_in(jkpt)%kpg_k_gather,1) 1134 sz2=size(bandfft_kpt_in(jkpt)%kpg_k_gather,2) 1135 ABI_MALLOC(bandfft_kpt_out(jkpt)%kpg_k_gather,(sz1,sz2)) 1136 bandfft_kpt_out(jkpt)%kpg_k_gather= & 1137 & bandfft_kpt_in(jkpt)%kpg_k_gather 1138 end if 1139 bandfft_kpt_out(jkpt)%flag3_is_allocated=bandfft_kpt_in(jkpt)%flag3_is_allocated 1140 if (allocated(bandfft_kpt_in(jkpt)%kg_k_gather_sym)) then 1141 sz1=size(bandfft_kpt_in(jkpt)%kg_k_gather_sym,1) 1142 sz2=size(bandfft_kpt_in(jkpt)%kg_k_gather_sym,2) 1143 ABI_MALLOC(bandfft_kpt_out(jkpt)%kg_k_gather_sym,(sz1,sz2)) 1144 bandfft_kpt_out(jkpt)%kg_k_gather_sym= & 1145 & bandfft_kpt_in(jkpt)%kg_k_gather_sym 1146 end if 1147 if (allocated(bandfft_kpt_in(jkpt)%rdispls_sym)) then 1148 sz1=size(bandfft_kpt_in(jkpt)%rdispls_sym) 1149 ABI_MALLOC(bandfft_kpt_out(jkpt)%rdispls_sym,(sz1)) 1150 bandfft_kpt_out(jkpt)%rdispls_sym= & 1151 & bandfft_kpt_in(jkpt)%rdispls_sym 1152 end if 1153 if (allocated(bandfft_kpt_in(jkpt)%recvcounts_sym)) then 1154 sz1=size(bandfft_kpt_in(jkpt)%recvcounts_sym) 1155 ABI_MALLOC(bandfft_kpt_out(jkpt)%recvcounts_sym,(sz1)) 1156 bandfft_kpt_out(jkpt)%recvcounts_sym= & 1157 & bandfft_kpt_in(jkpt)%recvcounts_sym 1158 end if 1159 if (allocated(bandfft_kpt_in(jkpt)%recvcounts_sym_tot)) then 1160 sz1=size(bandfft_kpt_in(jkpt)%recvcounts_sym_tot) 1161 ABI_MALLOC(bandfft_kpt_out(jkpt)%recvcounts_sym_tot,(sz1)) 1162 bandfft_kpt_out(jkpt)%recvcounts_sym_tot= & 1163 & bandfft_kpt_in(jkpt)%recvcounts_sym_tot 1164 end if 1165 if (allocated(bandfft_kpt_in(jkpt)%sdispls_sym)) then 1166 sz1=size(bandfft_kpt_in(jkpt)%sdispls_sym) 1167 ABI_MALLOC(bandfft_kpt_out(jkpt)%sdispls_sym,(sz1)) 1168 bandfft_kpt_out(jkpt)%sdispls_sym= & 1169 & bandfft_kpt_in(jkpt)%sdispls_sym 1170 end if 1171 if (allocated(bandfft_kpt_in(jkpt)%sendcounts_sym)) then 1172 sz1=size(bandfft_kpt_in(jkpt)%sendcounts_sym) 1173 ABI_MALLOC(bandfft_kpt_out(jkpt)%sendcounts_sym,(sz1)) 1174 bandfft_kpt_out(jkpt)%sendcounts_sym= & 1175 & bandfft_kpt_in(jkpt)%sendcounts_sym 1176 end if 1177 if (allocated(bandfft_kpt_in(jkpt)%sendcounts_sym_all)) then 1178 sz1=size(bandfft_kpt_in(jkpt)%sendcounts_sym_all) 1179 ABI_MALLOC(bandfft_kpt_out(jkpt)%sendcounts_sym_all,(sz1)) 1180 bandfft_kpt_out(jkpt)%sendcounts_sym_all= & 1181 & bandfft_kpt_in(jkpt)%sendcounts_sym_all 1182 end if 1183 if (allocated(bandfft_kpt_in(jkpt)%tab_proc)) then 1184 sz1=size(bandfft_kpt_in(jkpt)%tab_proc) 1185 ABI_MALLOC(bandfft_kpt_out(jkpt)%tab_proc,(sz1)) 1186 bandfft_kpt_out(jkpt)%tab_proc= & 1187 bandfft_kpt_in(jkpt)%tab_proc 1188 end if 1189 if (bandfft_kpt_in(jkpt)%have_to_reequilibrate) then 1190 bandfft_kpt_out(jkpt)%have_to_reequilibrate = .true. 1191 bandfft_kpt_out(jkpt)%npw_fft = bandfft_kpt_in(jkpt)%npw_fft 1192 sz1=size(bandfft_kpt_in(jkpt)%indices_pw_fft) 1193 ABI_MALLOC(bandfft_kpt_out(jkpt)%indices_pw_fft,(sz1)) 1194 bandfft_kpt_out(jkpt)%indices_pw_fft=bandfft_kpt_in(jkpt)%indices_pw_fft 1195 sz1=size(bandfft_kpt_in(jkpt)%sendcount_fft) 1196 ABI_MALLOC(bandfft_kpt_out(jkpt)%sendcount_fft,(sz1)) 1197 bandfft_kpt_out(jkpt)%sendcount_fft=bandfft_kpt_in(jkpt)%sendcount_fft 1198 sz1=size(bandfft_kpt_in(jkpt)%senddisp_fft) 1199 ABI_MALLOC(bandfft_kpt_out(jkpt)%senddisp_fft,(sz1)) 1200 bandfft_kpt_out(jkpt)%senddisp_fft=bandfft_kpt_in(jkpt)%senddisp_fft 1201 sz1=size(bandfft_kpt_in(jkpt)%recvcount_fft) 1202 ABI_MALLOC(bandfft_kpt_out(jkpt)%recvcount_fft,(sz1)) 1203 bandfft_kpt_out(jkpt)%recvcount_fft=bandfft_kpt_in(jkpt)%recvcount_fft 1204 sz1=size(bandfft_kpt_in(jkpt)%recvdisp_fft) 1205 ABI_MALLOC(bandfft_kpt_out(jkpt)%recvdisp_fft,(sz1)) 1206 bandfft_kpt_out(jkpt)%recvdisp_fft=bandfft_kpt_in(jkpt)%recvdisp_fft 1207 sz1=size(bandfft_kpt_in(jkpt)%kg_k_fft,1) 1208 sz2=size(bandfft_kpt_in(jkpt)%kg_k_fft,2) 1209 ABI_MALLOC(bandfft_kpt_out(jkpt)%kg_k_fft,(sz1,sz2)) 1210 bandfft_kpt_out(jkpt)%kg_k_fft=bandfft_kpt_in(jkpt)%kg_k_fft 1211 end if 1212 end do 1213 end do 1214 end if 1215 end if 1216 1217 end subroutine bandfft_kpt_copy
m_bandfft_kpt/bandfft_kpt_destroy [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_destroy
FUNCTION
Destroy a bandfft_kpt datastructure
INPUTS
OUTPUT
SIDE EFFECTS
bandfft_kpt_in=the datastructure to destroy
SOURCE
808 subroutine bandfft_kpt_destroy(bandfft_kpt_in) 809 810 !Arguments ------------------------------------ 811 type(bandfft_kpt_type) :: bandfft_kpt_in 812 !Local variables------------------------------- 813 814 ! *********************************************************************** 815 816 bandfft_kpt_in%flag1_is_allocated=0 817 bandfft_kpt_in%flag2_is_allocated=0 818 bandfft_kpt_in%flag3_is_allocated=0 819 bandfft_kpt_in%have_to_reequilibrate=.false. 820 821 #if defined HAVE_GPU && defined HAVE_YAKL 822 if (associated(bandfft_kpt_in%kg_k_gather)) then 823 if(bandfft_kpt_in%gpu_option==ABI_GPU_KOKKOS) then 824 ABI_FREE_MANAGED(bandfft_kpt_in%kg_k_gather) 825 else 826 ABI_FREE(bandfft_kpt_in%kg_k_gather) 827 end if 828 end if 829 #else 830 if (allocated(bandfft_kpt_in%kg_k_gather)) then 831 ABI_FREE(bandfft_kpt_in%kg_k_gather) 832 end if 833 #endif 834 835 if (allocated(bandfft_kpt_in%gbound)) then 836 ABI_FREE(bandfft_kpt_in%gbound) 837 end if 838 if (allocated(bandfft_kpt_in%recvcounts)) then 839 ABI_FREE(bandfft_kpt_in%recvcounts) 840 end if 841 if (allocated(bandfft_kpt_in%sendcounts)) then 842 ABI_FREE(bandfft_kpt_in%sendcounts) 843 end if 844 if (allocated(bandfft_kpt_in%rdispls)) then 845 ABI_FREE(bandfft_kpt_in%rdispls) 846 end if 847 if (allocated(bandfft_kpt_in%sdispls)) then 848 ABI_FREE(bandfft_kpt_in%sdispls) 849 end if 850 if (allocated(bandfft_kpt_in%ffnl_gather)) then 851 ABI_FREE(bandfft_kpt_in%ffnl_gather) 852 end if 853 854 #if defined HAVE_GPU && defined HAVE_YAKL 855 if (associated(bandfft_kpt_in%kinpw_gather)) then 856 if(bandfft_kpt_in%gpu_option==ABI_GPU_KOKKOS) then 857 ABI_FREE_MANAGED(bandfft_kpt_in%kinpw_gather) 858 else 859 ABI_FREE(bandfft_kpt_in%kinpw_gather) 860 end if 861 end if 862 #else 863 if (allocated(bandfft_kpt_in%kinpw_gather)) then 864 ABI_FREE(bandfft_kpt_in%kinpw_gather) 865 end if 866 #endif 867 868 if (allocated(bandfft_kpt_in%kpg_k_gather)) then 869 ABI_FREE(bandfft_kpt_in%kpg_k_gather) 870 end if 871 if (allocated(bandfft_kpt_in%ph3d_gather)) then 872 ABI_FREE(bandfft_kpt_in%ph3d_gather) 873 end if 874 if (allocated(bandfft_kpt_in%kg_k_gather_sym)) then 875 ABI_FREE(bandfft_kpt_in%kg_k_gather_sym) 876 end if 877 if (allocated(bandfft_kpt_in%rdispls_sym)) then 878 ABI_FREE(bandfft_kpt_in%rdispls_sym) 879 end if 880 if (allocated(bandfft_kpt_in%recvcounts_sym)) then 881 ABI_FREE(bandfft_kpt_in%recvcounts_sym) 882 end if 883 if (allocated(bandfft_kpt_in%recvcounts_sym_tot)) then 884 ABI_FREE(bandfft_kpt_in%recvcounts_sym_tot) 885 end if 886 if (allocated(bandfft_kpt_in%sdispls_sym)) then 887 ABI_FREE(bandfft_kpt_in%sdispls_sym) 888 end if 889 if (allocated(bandfft_kpt_in%sendcounts_sym)) then 890 ABI_FREE(bandfft_kpt_in%sendcounts_sym) 891 end if 892 if (allocated(bandfft_kpt_in%sendcounts_sym_all)) then 893 ABI_FREE(bandfft_kpt_in%sendcounts_sym_all) 894 end if 895 if (allocated(bandfft_kpt_in%tab_proc)) then 896 ABI_FREE(bandfft_kpt_in%tab_proc) 897 end if 898 if (allocated(bandfft_kpt_in%indices_pw_fft)) then 899 ABI_FREE(bandfft_kpt_in%indices_pw_fft) 900 end if 901 if (allocated(bandfft_kpt_in%sendcount_fft)) then 902 ABI_FREE(bandfft_kpt_in%sendcount_fft) 903 end if 904 if (allocated(bandfft_kpt_in%senddisp_fft)) then 905 ABI_FREE(bandfft_kpt_in%senddisp_fft) 906 end if 907 if (allocated(bandfft_kpt_in%recvcount_fft)) then 908 ABI_FREE(bandfft_kpt_in%recvcount_fft) 909 end if 910 if (allocated(bandfft_kpt_in%recvdisp_fft)) then 911 ABI_FREE(bandfft_kpt_in%recvdisp_fft) 912 end if 913 if (allocated(bandfft_kpt_in%kg_k_fft)) then 914 ABI_FREE(bandfft_kpt_in%kg_k_fft) 915 end if 916 917 end subroutine bandfft_kpt_destroy
m_bandfft_kpt/bandfft_kpt_destroy_array [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_destroy_array
FUNCTION
Clean and destroy an array of bandfft_kpt datastructures
INPUTS
OUTPUT
SIDE EFFECTS
bandfft_kpt_in(:)=the array of datastructure to destroy
SOURCE
939 subroutine bandfft_kpt_destroy_array(bandfft_kpt_in,mpi_enreg) 940 941 !Arguments ------------------------------------ 942 type(bandfft_kpt_type),pointer :: bandfft_kpt_in(:) 943 type(MPI_type), intent(inout) :: mpi_enreg 944 945 !Local variables------------------------------- 946 integer :: ikpt_this_proc,isppol,ikpt,mkmem,nband,nkpt,nsppol 947 character(len=500) :: msg 948 949 ! *********************************************************************** 950 951 if (xmpi_paral==0) return 952 953 if (associated(bandfft_kpt_in)) then 954 mkmem =size(bandfft_kpt_in) 955 nkpt=size(mpi_enreg%proc_distrb,1) 956 nband=size(mpi_enreg%proc_distrb,2) 957 nsppol=size(mpi_enreg%proc_distrb,3) 958 if (nsppol==0.or.nkpt==0) then 959 msg=' mpi_enreg%proc_distrb should be allocated !' 960 ABI_BUG(msg) 961 end if 962 nkpt=size(mpi_enreg%my_kpttab) 963 if (nkpt==0) then 964 msg=' mpi_enreg%my_kpttab should be allocated !' 965 ABI_BUG(msg) 966 end if 967 do isppol=1,nsppol 968 do ikpt=1,nkpt 969 if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband,isppol,mpi_enreg%me_kpt)) cycle 970 ikpt_this_proc=mpi_enreg%my_kpttab(ikpt) 971 if ((ikpt_this_proc>mkmem) .or.(ikpt_this_proc<=0)) then 972 msg=' The bandfft tab cannot be deallocated !' 973 ABI_BUG(msg) 974 end if 975 call bandfft_kpt_destroy(bandfft_kpt_in(ikpt_this_proc)) 976 end do 977 end do 978 ABI_FREE(bandfft_kpt_in) 979 nullify(bandfft_kpt_in) 980 end if 981 982 end subroutine bandfft_kpt_destroy_array
m_bandfft_kpt/bandfft_kpt_get_ikpt [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_get_ikpt
FUNCTION
Get the value of bandfft_kpt_current_ikpt variable, i.e. current index of bandfft_kpt loaded in memory INPUT
OUTPUT
bandfft_kpt_get_ikpt= current index of bandfft_kpt
SOURCE
2168 function bandfft_kpt_get_ikpt() 2169 2170 !Arguments ------------------------------- 2171 integer :: bandfft_kpt_get_ikpt 2172 !Local variables------------------------------- 2173 2174 ! ********************************************************************* 2175 2176 bandfft_kpt_get_ikpt=bandfft_kpt_current_ikpt 2177 2178 end function bandfft_kpt_get_ikpt
m_bandfft_kpt/bandfft_kpt_init1 [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_init1
FUNCTION
Init all (or part of) scalars and pointers in a bandfft_kpt datastructure
INPUTS
istwfk(nkpt) = input option parameter that describes the storage of wfs kg(3,mpw*mkmem) = dimensionless coords of G vecs in basis sphere at k point mgfft = maximum single fft dimension (IN) mkmem = number of k points which can fit in memory; set to 0 if use disk mpi_enreg = information about MPI parallelization mpw = maximum number of planewaves as dimensioned in calling routine nband(nkpt*nsppol) = number of bands at each k point nkpt = number of k points npwarr(nkpt) = array holding npw for each k point, taking into account the effect of istwfk, and the spreading over processors nsppol = 1 for unpolarized, 2 for polarized
OUTPUT
within the bandfft_kpt data_type : Initialize and compute gbound = sphere boundary info idatarecv0 = position of the planewave coordinates (0,0,0) istwf_k = input option parameter that describes the storage of wfs kg_k_gather = planewave coordinates (of the processor + sended by other processors band) kg_k_gather_sym = planewave coordinates (kg_k_gather + opposited planewave coordinates sended by the processors fft) ndatarecv = total number of values received by the processor and sended by the other processors band ndatasend_sym = number of sended values to the processors fft to create opposited planewave coordinates ndatarecv_tot = total number of received values by the processor (ndatarecv + number of received opposited planewave coordinates) recvcounts = number of values received by the processor from each processor band recvcounts_sym = number of values received by the processor from each processor fft recvcounts_sym_tot = number of values received by each processor from the other processors fft rdispls = positions of values received by the processor from each processor band rdispls_sym = positions of values received by the processor from each processor fft sendcounts = number of values sended by the processor to each processor band sendcounts_sym = number of values sended by the processor to each processor fft sendcounts_sym_all = number of values sended by each processor to the other processors fft sdispls = postions of values sended by the processor to each processor band sdispls_sym = postions of values sended by the processor to each processor fft tab_proc = positions of opposited planewave coordinates in the list of the processors fft
SIDE EFFECTS
bandfft_kpt_in=<type(bandfft_kpt)>=bandfft_kpt datastructure
SOURCE
213 subroutine bandfft_kpt_init1(bandfft_kpt_in,istwfk,kg,mgfft,mkmem,mpi_enreg,mpw,nband,nkpt,npwarr,nsppol,gpu_option) 214 215 !Arguments ------------------------------------ 216 !scalars 217 integer,intent(in) :: mgfft,mkmem,mpw,nkpt,nsppol 218 type(bandfft_kpt_type),pointer :: bandfft_kpt_in(:) 219 type(MPI_type),intent(inout) :: mpi_enreg 220 integer,intent(in),optional :: gpu_option 221 !arrays 222 integer,intent(in) :: istwfk(nkpt),nband(nkpt*nsppol) 223 integer,intent(in) :: kg(3,mpw*mkmem),npwarr(nkpt) 224 225 !Local variables------------------------------- 226 !scalars 227 integer :: comm_band,comm_fft,idatarecv,idatarecv0,ierr,ikg,ikpt,ikpt_this_proc,iproc,isppol,istwf_k,itest,jsendloc 228 integer :: me_fft,me_kpt,n2,nband_k,ndatarecv,ndatarecv_tot,ndatasend_sym,nproc_band,nproc_fft,npw_fft,npw_k,npw_tot 229 logical :: reequilibrate_allocated 230 character(len=500) :: message 231 !arrays 232 integer,allocatable :: buff_kg(:,:),gbound(:,:),indices_pw_fft(:),kg_k(:,:),kg_k_fft(:,:),kg_k_gather(:,:) 233 integer,allocatable :: kg_k_gather_all(:,:),kg_k_gather_send(:,:),kg_k_gather_sym(:,:) 234 integer,allocatable :: npw_per_proc(:) 235 integer,allocatable :: rdispls(:),rdispls_all(:),rdispls_sym(:),rdispls_sym_loc(:) 236 integer,allocatable :: recvcounts(:),recvcounts_sym(:),recvcounts_fft(:),recvcounts_sym_loc(:),recvcounts_sym_tot(:) 237 integer,allocatable :: recvdisp_fft(:),sdispls(:),sdispls_sym_loc(:),sdispls_sym(:) 238 integer,allocatable :: sendcounts(:),sendcounts_fft(:),sendcounts_sym(:),sendcounts_sym_all(:),sendcounts_sym_loc(:) 239 integer,allocatable :: senddisp_fft(:),sum_kg(:),tab_proc(:) 240 241 ! ********************************************************************* 242 243 If(mpi_enreg%paral_kgb/=1) then 244 nullify(bandfft_kpt_in) 245 return 246 end if 247 248 !--------------------------------------------- 249 !Initialisation 250 !--------------------------------------------- 251 nproc_fft = mpi_enreg%nproc_fft 252 nproc_band = mpi_enreg%nproc_band 253 254 me_fft = mpi_enreg%me_fft 255 me_kpt = mpi_enreg%me_kpt 256 257 comm_band = mpi_enreg%comm_band 258 comm_fft = mpi_enreg%comm_fft 259 260 !============================================================================= 261 !Compute and store various tabs in bandfft_kpt(ikpt) data_struc 262 !These ones will be used in following subroutines: 263 !vtorho, mkrho, prep_nonlop, prep_fourwf, prep_getghc... 264 !============================================================================= 265 266 ABI_MALLOC(bandfft_kpt_in,(mkmem)) 267 268 ABI_MALLOC(sdispls ,(nproc_band)) 269 ABI_MALLOC(sendcounts ,(nproc_band)) 270 ABI_MALLOC(rdispls ,(nproc_band)) 271 ABI_MALLOC(recvcounts ,(nproc_band)) 272 ABI_MALLOC(sendcounts_fft,(nproc_fft)) 273 ABI_MALLOC(senddisp_fft ,(nproc_fft)) 274 ABI_MALLOC(recvcounts_fft,(nproc_fft)) 275 ABI_MALLOC(recvdisp_fft ,(nproc_fft)) 276 277 bandfft_kpt_in(:)%flag1_is_allocated=0 278 bandfft_kpt_in(:)%flag2_is_allocated=0 279 bandfft_kpt_in(:)%flag3_is_allocated=0 280 bandfft_kpt_in(:)%gpu_option=ABI_GPU_DISABLED 281 if(present(gpu_option)) bandfft_kpt_in(:)%gpu_option=gpu_option 282 283 do isppol=1,nsppol 284 ikg=0 285 do ikpt=1,nkpt 286 npw_k=npwarr(ikpt) 287 istwf_k=istwfk(ikpt) 288 nband_k=nband(ikpt+(isppol-1)*nkpt) 289 if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,isppol,me_kpt))cycle 290 ikpt_this_proc=mpi_enreg%my_kpttab(ikpt) 291 if ((ikpt_this_proc > mkmem).or.(ikpt_this_proc==0)) then 292 message = ' this bandfft tab is not allocated !' 293 ABI_ERROR(message) 294 end if 295 if (bandfft_kpt_in(ikpt_this_proc)%flag1_is_allocated==0) then 296 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%gbound ,(2*mgfft+8,2)) 297 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%recvcounts,(nproc_band)) 298 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sendcounts,(nproc_band)) 299 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%rdispls ,(nproc_band)) 300 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sdispls ,(nproc_band)) 301 bandfft_kpt_in(ikpt_this_proc)%flag1_is_allocated=1 302 end if 303 304 ! Initialize various quantities for the reequilibration step 305 reequilibrate_allocated=(isppol==2.and.bandfft_kpt_in(ikpt_this_proc)%have_to_reequilibrate) 306 bandfft_kpt_in(ikpt_this_proc)%have_to_reequilibrate = .false. 307 bandfft_kpt_in(ikpt_this_proc)%npw_fft = 0 308 309 call xmpi_allgather(npw_k,recvcounts,comm_band,ierr) 310 rdispls(1)=0 311 do iproc=2,nproc_band 312 rdispls(iproc)=rdispls(iproc-1)+recvcounts(iproc-1) 313 end do 314 ndatarecv=rdispls(nproc_band)+recvcounts(nproc_band) 315 316 ABI_MALLOC(kg_k_gather,(3,ndatarecv)) 317 ABI_MALLOC(kg_k,(3,npw_k)) 318 kg_k(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg) 319 call xmpi_allgatherv(kg_k,3*npw_k,kg_k_gather,3*recvcounts(:),3*rdispls(:),comm_band,ierr) 320 321 sendcounts(:)=npw_k*mpi_enreg%bandpp 322 do iproc=1,nproc_band 323 sdispls(iproc)=(iproc-1)*npw_k*mpi_enreg%bandpp 324 end do 325 326 ! ============================================================================ 327 ! Here we compute gbound, as well for istwf_k=1 as for istwf_k=2 and store it 328 ! ============================================================================ 329 !MG: Why don't we compute gbound locally by just computing the full G-sphere from ecut 330 ABI_MALLOC(npw_per_proc,(nproc_fft)) 331 ABI_MALLOC(rdispls_all,(nproc_fft)) 332 ABI_MALLOC(gbound,(2*mgfft+8,2)) 333 if (mgfft>0) gbound(:,:)=0 334 if (istwf_k==1) then 335 call xmpi_allgather(ndatarecv,npw_per_proc,mpi_enreg%comm_fft,ierr) 336 rdispls_all(1)=0 337 do iproc=2,nproc_fft 338 rdispls_all(iproc)=rdispls_all(iproc-1)+npw_per_proc(iproc-1) 339 end do 340 npw_tot=rdispls_all(nproc_fft)+npw_per_proc(nproc_fft) 341 ABI_MALLOC(kg_k_gather_all,(3,npw_tot)) 342 call xmpi_allgatherv(kg_k_gather,& 343 & 3*ndatarecv,kg_k_gather_all,3*npw_per_proc(:),3*rdispls_all,mpi_enreg%comm_fft,ierr) 344 if (mgfft>0) then 345 call sphereboundary(gbound,istwf_k,kg_k_gather_all,mgfft,npw_tot) 346 end if 347 348 else if (istwf_k==2) then 349 350 ! ============================================================================ 351 ! In this case, we have to add the opposite values in the kg_k_gather tab 352 ! before computing gbound 353 ! ============================================================================ 354 355 ! Allocation 356 ABI_MALLOC(tab_proc ,(ndatarecv)) 357 ABI_MALLOC(sendcounts_sym ,(nproc_fft)) 358 ABI_MALLOC(sendcounts_sym_all,(nproc_fft*nproc_fft)) 359 ABI_MALLOC(sdispls_sym ,(nproc_fft)) 360 ABI_MALLOC(recvcounts_sym ,(nproc_fft)) 361 ABI_MALLOC(recvcounts_sym_tot,(nproc_fft)) 362 ABI_MALLOC(rdispls_sym ,(nproc_fft)) 363 ABI_MALLOC(sendcounts_sym_loc,(nproc_fft)) 364 ABI_MALLOC(sdispls_sym_loc ,(nproc_fft)) 365 ABI_MALLOC(recvcounts_sym_loc,(nproc_fft)) 366 ABI_MALLOC(rdispls_sym_loc ,(nproc_fft)) 367 368 ! Initialisation 369 tab_proc(:) = 0 370 sendcounts_sym(:) = 0 371 sendcounts_sym_all(:) = 0 372 sdispls_sym(:) = 0 373 recvcounts_sym(:) = 0 374 recvcounts_sym_tot(:) = 0 375 376 ! Localisation of kg_k==[0 0 0] 377 ABI_MALLOC(sum_kg,(ndatarecv)) 378 idatarecv0 = -1 379 ndatasend_sym = ndatarecv 380 sum_kg=sum(abs(kg_k_gather),1) 381 if (count(sum_kg==0)/=0) then 382 do idatarecv=1,ndatarecv 383 if (sum_kg(idatarecv)==0) idatarecv0=idatarecv 384 end do 385 ndatasend_sym = ndatarecv-1 386 end if 387 388 ! Localisation of the processor where the vector -k2 is 389 n2 = mpi_enreg%distribfft%n2_coarse 390 do idatarecv=1,ndatarecv 391 if (idatarecv/=idatarecv0) then 392 tab_proc(idatarecv) = mpi_enreg%distribfft%tab_fftwf2_distrib(modulo(-kg_k_gather(2,idatarecv),n2) + 1) 393 else 394 tab_proc(idatarecv) = -1 395 end if 396 end do 397 398 ! Number of values send by the processor to the others 399 do iproc=1,nproc_fft 400 sendcounts_sym(iproc) = count(tab_proc(:)==(iproc-1)) 401 end do 402 403 ! Save sendcounts_sym for each processor in sendcounts_sym_all 404 ! knowed by all processors of comm_fft 405 rdispls_sym(1)=0 406 do iproc=2,nproc_fft 407 rdispls_sym(iproc)= nproc_fft*(iproc-1) 408 end do 409 recvcounts_sym(:)=nproc_fft 410 call xmpi_allgatherv(sendcounts_sym(:),nproc_fft,& 411 & sendcounts_sym_all(:),recvcounts_sym,rdispls_sym,comm_fft,ierr) 412 413 ! Calculation of the dimension of kg_k_gather_sym for each processor 414 ! recvcounts_sym_tot is knowed by all processors of comm_fft 415 call xmpi_sum(sendcounts_sym,recvcounts_sym_tot,nproc_fft,comm_fft,ierr) 416 417 ! Dimension of kg_k_gather_sym 418 ndatarecv_tot = ndatarecv+recvcounts_sym_tot(me_fft+1) 419 420 ! Intialize kg_k_gather_sym 421 ABI_MALLOC(kg_k_gather_sym,(3,ndatarecv_tot)) 422 kg_k_gather_sym(:,:)=0 423 kg_k_gather_sym(:,1:ndatarecv) = kg_k_gather(:,:) 424 425 ! Allocation and initialisation 426 ABI_MALLOC(kg_k_gather_send,(3,ndatasend_sym)) 427 kg_k_gather_send(:,:)=0 428 429 ! The values are sorted in blocks 430 jsendloc=0 431 do iproc=1,nproc_fft 432 433 ! Position of the beginning of the block 434 sdispls_sym(iproc)=jsendloc 435 436 ! Creation of the blocks 437 do idatarecv=1,ndatarecv 438 if (tab_proc(idatarecv)==(iproc-1)) then 439 jsendloc=jsendloc+1 440 kg_k_gather_send(:,jsendloc) = -kg_k_gather(:,idatarecv) 441 end if 442 end do 443 end do 444 445 ! Position of received data 446 rdispls_sym(1)= ndatarecv 447 recvcounts_sym(1)= sendcounts_sym_all((me_fft+1)) 448 do iproc=2,nproc_fft 449 rdispls_sym(iproc) = rdispls_sym(iproc-1) + & 450 sendcounts_sym_all((me_fft+1)+(iproc-2)*nproc_fft) 451 recvcounts_sym(iproc) = sendcounts_sym_all((me_fft+1)+(iproc-1)*nproc_fft) 452 end do 453 454 ! Exchange of kg_k 455 sendcounts_sym_loc = sendcounts_sym*3 456 sdispls_sym_loc = sdispls_sym *3 457 recvcounts_sym_loc = recvcounts_sym*3 458 rdispls_sym_loc = rdispls_sym *3 459 call xmpi_alltoallv(kg_k_gather_send(:,:),sendcounts_sym_loc,sdispls_sym_loc,& 460 & kg_k_gather_sym(:,:) ,recvcounts_sym_loc,rdispls_sym_loc,comm_fft,ierr) 461 462 ! Store the following data in the bandfft_kpt_in data_struc 463 ikpt_this_proc=mpi_enreg%my_kpttab(ikpt) 464 if (bandfft_kpt_in(ikpt_this_proc)%flag3_is_allocated==0) then 465 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather_sym,(3,ndatarecv_tot)) 466 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%rdispls_sym,(nproc_fft)) 467 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%recvcounts_sym,(nproc_fft)) 468 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%recvcounts_sym_tot,(nproc_fft)) 469 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sdispls_sym,(nproc_fft)) 470 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sendcounts_sym,(nproc_fft)) 471 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sendcounts_sym_all,(nproc_fft*nproc_fft)) 472 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%tab_proc,(ndatarecv)) 473 bandfft_kpt_in(ikpt_this_proc)%flag3_is_allocated=1 474 end if 475 476 bandfft_kpt_in(ikpt_this_proc)%idatarecv0 =idatarecv0 477 bandfft_kpt_in(ikpt_this_proc)%ndatarecv_tot =ndatarecv_tot 478 bandfft_kpt_in(ikpt_this_proc)%ndatasend_sym =ndatasend_sym 479 bandfft_kpt_in(ikpt_this_proc)%kg_k_gather_sym(:,:) =kg_k_gather_sym(:,:) 480 bandfft_kpt_in(ikpt_this_proc)%rdispls_sym(:) =rdispls_sym(:) 481 bandfft_kpt_in(ikpt_this_proc)%recvcounts_sym(:) =recvcounts_sym(:) 482 bandfft_kpt_in(ikpt_this_proc)%recvcounts_sym_tot(:)=recvcounts_sym_tot(:) 483 bandfft_kpt_in(ikpt_this_proc)%sdispls_sym(:) =sdispls_sym(:) 484 bandfft_kpt_in(ikpt_this_proc)%sendcounts_sym(:) =sendcounts_sym(:) 485 bandfft_kpt_in(ikpt_this_proc)%sendcounts_sym_all(:)=sendcounts_sym_all(:) 486 bandfft_kpt_in(ikpt_this_proc)%tab_proc(:) =tab_proc(:) 487 488 ABI_FREE(tab_proc) 489 ABI_FREE(sendcounts_sym) 490 ABI_FREE(sendcounts_sym_all) 491 ABI_FREE(sdispls_sym) 492 ABI_FREE(recvcounts_sym) 493 ABI_FREE(recvcounts_sym_tot) 494 ABI_FREE(rdispls_sym) 495 ABI_FREE(kg_k_gather_sym) 496 ABI_FREE(sendcounts_sym_loc) 497 ABI_FREE(recvcounts_sym_loc) 498 ABI_FREE(sdispls_sym_loc) 499 ABI_FREE(rdispls_sym_loc) 500 ABI_FREE(kg_k_gather_send) 501 ABI_FREE(sum_kg) 502 503 ! Then compute gbound 504 call xmpi_allgather(ndatarecv_tot,npw_per_proc,mpi_enreg%comm_fft,ierr) 505 rdispls_all(1)=0 506 do iproc=2,nproc_fft 507 rdispls_all(iproc)=rdispls_all(iproc-1)+npw_per_proc(iproc-1) 508 end do 509 npw_tot=rdispls_all(nproc_fft)+npw_per_proc(nproc_fft) 510 ABI_MALLOC(kg_k_gather_all,(3,npw_tot)) 511 call xmpi_allgatherv(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather_sym,& 512 & 3*ndatarecv_tot,kg_k_gather_all,3*npw_per_proc(:),3*rdispls_all,mpi_enreg%comm_fft,ierr) 513 if (mgfft>0) then 514 call sphereboundary(gbound,istwf_k,kg_k_gather_all,mgfft,npw_tot) 515 end if 516 517 ! Only calculations with istwfk=1 or 2 518 else 519 write(message, '(a,i0,a)' )' the value istwfk=',istwf_k,' is not allowed in case of bandfft parallelization!' 520 ABI_BUG(message) 521 end if 522 ABI_FREE(kg_k_gather_all) 523 ABI_FREE(npw_per_proc) 524 ABI_FREE(rdispls_all) 525 ! ============================================================================ 526 ! End of gbound 527 ! ============================================================================ 528 529 ! Check if there is k_G vectors have been redistributed (after unbalancing dectecion) 530 ! Note that FFT load balancing is directly related to unbalancing detection 531 ! made in kpgsph routine. 532 itest=0 ; n2 = mpi_enreg%distribfft%n2_coarse 533 do idatarecv=1,ndatarecv 534 iproc=mpi_enreg%distribfft%tab_fftwf2_distrib(modulo(kg_k_gather(2,idatarecv),n2)+1) 535 if (iproc/=me_fft) itest=itest+1 536 end do 537 call xmpi_sum(itest,mpi_enreg%comm_fft,ierr) 538 if (itest>0) then 539 write(message, '(a,i4,3a)' ) & 540 & 'There is a load unbalancing for the FFT parallelization (kpt',ikpt,').',ch10,& 541 & 'Plane-wave components will be redistributed before each FFT!' 542 ABI_COMMENT(message) 543 end if 544 bandfft_kpt_in(ikpt_this_proc)%have_to_reequilibrate=(itest>0) 545 546 ! If yes, store relevant data 547 if(bandfft_kpt_in(ikpt_this_proc)%have_to_reequilibrate) then 548 n2 = mpi_enreg%distribfft%n2_coarse 549 sendcounts_fft(:) = 0 550 recvcounts_fft(:) = 0 551 senddisp_fft(:) = 0 552 recvdisp_fft(:) = 0 553 do idatarecv = 1 ,ndatarecv 554 iproc = mpi_enreg%distribfft%tab_fftwf2_distrib(modulo(kg_k_gather(2,idatarecv),n2) + 1) 555 sendcounts_fft(iproc + 1) = sendcounts_fft(iproc + 1) + 1 556 end do 557 call xmpi_alltoall(sendcounts_fft,1,recvcounts_fft,1,mpi_enreg%comm_fft,ierr) 558 do iproc =2, nproc_fft 559 senddisp_fft(iproc) = senddisp_fft(iproc - 1) + sendcounts_fft(iproc-1) 560 recvdisp_fft(iproc) = recvdisp_fft(iproc - 1) + recvcounts_fft(iproc-1) 561 end do 562 npw_fft = recvdisp_fft(nproc_fft) + recvcounts_fft(nproc_fft) ! nb plane wave for fourwf call 563 ABI_MALLOC(buff_kg,(3,ndatarecv)) ! for sorting kg_k 564 ABI_MALLOC(kg_k_fft,(3,npw_fft)) 565 ABI_MALLOC(indices_pw_fft,(ndatarecv)) 566 !filling of sorted send buffers 567 sendcounts_fft(:) = 0 568 do idatarecv = 1 ,ndatarecv 569 iproc = mpi_enreg%distribfft%tab_fftwf2_distrib(modulo(kg_k_gather(2,idatarecv),n2) + 1) 570 sendcounts_fft(iproc + 1) = sendcounts_fft(iproc + 1) + 1 571 indices_pw_fft(idatarecv) = senddisp_fft(iproc+1) + sendcounts_fft(iproc+1) 572 buff_kg(1:3,indices_pw_fft(idatarecv)) = kg_k_gather(1:3,idatarecv) 573 end do 574 575 call xmpi_alltoallv(buff_kg, 3*sendcounts_fft, 3*senddisp_fft, & 576 & kg_k_fft,3*recvcounts_fft, 3*recvdisp_fft, mpi_enreg%comm_fft,ierr) 577 578 if (.not.reequilibrate_allocated) then 579 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kg_k_fft, (3,npw_fft) ) 580 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%indices_pw_fft, (ndatarecv)) 581 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%sendcount_fft, (nproc_fft)) 582 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%senddisp_fft, (nproc_fft)) 583 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%recvcount_fft, (nproc_fft)) 584 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%recvdisp_fft, (nproc_fft)) 585 end if 586 bandfft_kpt_in(ikpt_this_proc)%npw_fft = npw_fft 587 bandfft_kpt_in(ikpt_this_proc)%indices_pw_fft(:) = indices_pw_fft(:) 588 bandfft_kpt_in(ikpt_this_proc)%sendcount_fft (:) = sendcounts_fft (:) 589 bandfft_kpt_in(ikpt_this_proc)%senddisp_fft(:) = senddisp_fft(:) 590 bandfft_kpt_in(ikpt_this_proc)%recvcount_fft(:) = recvcounts_fft(:) 591 bandfft_kpt_in(ikpt_this_proc)%recvdisp_fft(:) = recvdisp_fft(:) 592 bandfft_kpt_in(ikpt_this_proc)%kg_k_fft(:,:) = kg_k_fft(:,:) 593 ABI_FREE(buff_kg ) 594 ABI_FREE(kg_k_fft) 595 ABI_FREE(indices_pw_fft) 596 end if 597 598 ! Tabs which are common to istwf_k=1 and 2 599 #if defined HAVE_GPU && defined HAVE_YAKL 600 if (.not. associated(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather)) then 601 if(bandfft_kpt_in(ikpt_this_proc)%gpu_option==ABI_GPU_KOKKOS) then 602 ABI_MALLOC_MANAGED(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather,(/3,ndatarecv/)) 603 else 604 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather,(3,ndatarecv)) 605 end if 606 end if 607 #else 608 if (.not. allocated(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather)) then 609 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kg_k_gather,(3,ndatarecv)) 610 end if 611 #endif 612 613 bandfft_kpt_in(ikpt_this_proc)%recvcounts(:) =recvcounts(:) 614 bandfft_kpt_in(ikpt_this_proc)%sendcounts(:) =sendcounts(:) 615 bandfft_kpt_in(ikpt_this_proc)%rdispls(:) =rdispls(:) 616 bandfft_kpt_in(ikpt_this_proc)%sdispls(:) =sdispls(:) 617 bandfft_kpt_in(ikpt_this_proc)%gbound(:,:) =gbound(:,:) 618 bandfft_kpt_in(ikpt_this_proc)%kg_k_gather(:,:)=kg_k_gather(:,:) 619 bandfft_kpt_in(ikpt_this_proc)%ndatarecv =ndatarecv 620 bandfft_kpt_in(ikpt_this_proc)%istwf_k =istwf_k 621 bandfft_kpt_in(ikpt_this_proc)%npw_tot =npw_tot 622 ABI_FREE(kg_k_gather) 623 ABI_FREE(kg_k) 624 ABI_FREE(gbound) 625 626 ikg=ikg+npw_k 627 end do 628 end do 629 ABI_FREE(recvcounts) 630 ABI_FREE(sendcounts) 631 ABI_FREE(rdispls) 632 ABI_FREE(sdispls) 633 ABI_FREE(sendcounts_fft) 634 ABI_FREE(senddisp_fft) 635 ABI_FREE(recvcounts_fft) 636 ABI_FREE(recvdisp_fft) 637 638 !============================================================================= 639 !End of computation and storage of the bandfft_kpt(ikpt) data_struc 640 !============================================================================= 641 642 end subroutine bandfft_kpt_init1
m_bandfft_kpt/bandfft_kpt_init2 [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_init2
FUNCTION
Init part of scalars and pointers in a bandfft_kpt datastructure
INPUTS
dimffnl=second dimension of ffnl (1+number of derivatives) ffnl_gather(ndatarecv,dimffnl,lmnmax,ntypat)=nonlocal form factors on basis sphere. kinpw_gather(:)=(modified) kinetic energy for each plane wave (Hartree) kpg_k_gather(ndatarecv,nkpg)=k+G vector for a given k point lmnmax=if useylm=1, max number of (l,m,n) comp. over all type of psps =if useylm=0, max number of (l,n) comp. over all type of psps matblk=dimension of the array ph3d mkmem =number of k points which can fit in memory; set to 0 if use disk ndatarecv= dimension of the arrays nkpg=second dimension of kpg_k (0 if useylm=0) ntypat=number of types of atoms in unit cell. ph3d_gather(2,ndatarecv,matblk)=3-dim structure factors, for each atom and plane wave.
OUTPUT
SIDE EFFECTS
SOURCE
674 subroutine bandfft_kpt_init2(bandfft_kpt_in,dimffnl,ffnl_gather,ikpt_this_proc,kinpw_gather,& 675 & kpg_k_gather,lmnmax,matblk,mkmem,ndatarecv,nkpg,ntypat,ph3d_gather) 676 677 !Arguments ------------------------------- 678 integer, intent(in) :: dimffnl,ikpt_this_proc,lmnmax,matblk,mkmem,ndatarecv,nkpg,ntypat 679 !Local variables------------------------------- 680 real(dp),intent(in) :: ffnl_gather(:,:,:,:),kinpw_gather(:) 681 real(dp),intent(in) :: kpg_k_gather(:,:),ph3d_gather(:,:,:) 682 type(bandfft_kpt_type), intent(inout) :: bandfft_kpt_in (mkmem) 683 684 ! ********************************************************************* 685 686 if (allocated(bandfft_kpt_in(ikpt_this_proc)%ffnl_gather)) then 687 ABI_FREE(bandfft_kpt_in(ikpt_this_proc)%ffnl_gather) 688 end if 689 if (size(ffnl_gather)>0) then 690 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%ffnl_gather,(ndatarecv,dimffnl,lmnmax,ntypat)) 691 bandfft_kpt_in(ikpt_this_proc)%ffnl_gather(:,:,:,:)=ffnl_gather(:,:,:,:) 692 else 693 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%ffnl_gather,(0,0,0,0)) 694 end if 695 696 if (allocated(bandfft_kpt_in(ikpt_this_proc)%ph3d_gather)) then 697 ABI_FREE(bandfft_kpt_in(ikpt_this_proc)%ph3d_gather) 698 end if 699 if (size(ph3d_gather)>0) then 700 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%ph3d_gather,(2,ndatarecv,matblk)) 701 bandfft_kpt_in(ikpt_this_proc)%ph3d_gather(:,:,:) =ph3d_gather(:,:,:) 702 else 703 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%ph3d_gather,(0,0,0)) 704 end if 705 706 if (allocated(bandfft_kpt_in(ikpt_this_proc)%kpg_k_gather)) then 707 ABI_FREE(bandfft_kpt_in(ikpt_this_proc)%kpg_k_gather) 708 end if 709 if (size(kpg_k_gather)>0) then 710 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kpg_k_gather,(ndatarecv,nkpg)) 711 bandfft_kpt_in(ikpt_this_proc)%kpg_k_gather(:,:) =kpg_k_gather(:,:) 712 else 713 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kpg_k_gather,(0,0)) 714 end if 715 716 #if defined HAVE_GPU && defined HAVE_YAKL 717 if (associated(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather)) then 718 if(bandfft_kpt_in(ikpt_this_proc)%gpu_option==ABI_GPU_KOKKOS) then 719 ABI_FREE_MANAGED(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather) 720 else 721 ABI_FREE(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather) 722 end if 723 end if 724 if(bandfft_kpt_in(ikpt_this_proc)%gpu_option==ABI_GPU_KOKKOS) then 725 if (size(kinpw_gather)>0) then 726 ABI_MALLOC_MANAGED(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(/ndatarecv/)) 727 bandfft_kpt_in(ikpt_this_proc)%kinpw_gather(:) =kinpw_gather(:) 728 else 729 ABI_MALLOC_MANAGED(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(/0/)) 730 end if 731 else 732 if (size(kinpw_gather)>0) then 733 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(ndatarecv)) 734 bandfft_kpt_in(ikpt_this_proc)%kinpw_gather(:) =kinpw_gather(:) 735 else 736 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(0)) 737 end if 738 end if 739 #else 740 if (allocated(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather)) then 741 ABI_FREE(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather) 742 end if 743 if (size(kinpw_gather)>0) then 744 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(ndatarecv)) 745 bandfft_kpt_in(ikpt_this_proc)%kinpw_gather(:) =kinpw_gather(:) 746 else 747 ABI_MALLOC(bandfft_kpt_in(ikpt_this_proc)%kinpw_gather,(0)) 748 end if 749 #endif 750 751 bandfft_kpt_in(ikpt_this_proc)%flag2_is_allocated=1 752 753 end subroutine bandfft_kpt_init2
m_bandfft_kpt/bandfft_kpt_reset [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_reset
FUNCTION
Reset flags a bandfft_kpt datastructure
INPUTS
OUTPUT
SIDE EFFECTS
bandfft_kpt_in=the datastructure to nullify
SOURCE
774 subroutine bandfft_kpt_reset(bandfft_kpt_in) 775 776 !Arguments ------------------------------------ 777 type(bandfft_kpt_type) :: bandfft_kpt_in 778 !Local variables------------------------------- 779 780 ! *********************************************************************** 781 782 bandfft_kpt_in%flag1_is_allocated=0 783 bandfft_kpt_in%flag2_is_allocated=0 784 bandfft_kpt_in%flag3_is_allocated=0 785 bandfft_kpt_in%have_to_reequilibrate=.false. 786 787 end subroutine bandfft_kpt_reset
m_bandfft_kpt/bandfft_kpt_restoretabs [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_restoretabs
FUNCTION
Restore some (kpt-dependent) tabs into a bandfft_kpt datastructure INPUT
OUTPUT
=== Arrays to be eventually restored ===== ffnl(:,:,:,:)=nonlocal form factors on basis sphere ph3d(:,:,:)=3-dim structure factors, for each atom and plane wave kpg(:,:)=k+G vector for a given k point kinpw(:)=kinetic energy for each plane wave (Hartree)
SIDE EFFECTS
bandfft_kpt_out=<type(bandfft_kpt)>=bandfft_kpt datastructure
SOURCE
2025 subroutine bandfft_kpt_restoretabs(bandfft_kpt_out,ffnl,ph3d,kpg,kinpw) 2026 2027 !Arguments ------------------------------- 2028 type(bandfft_kpt_type),intent(inout) :: bandfft_kpt_out 2029 real(dp),intent(inout),allocatable,optional :: ffnl(:,:,:,:),ph3d(:,:,:),kpg(:,:),kinpw(:) 2030 !Local variables------------------------------- 2031 integer :: is1,is2,is3,is4 2032 2033 ! ********************************************************************* 2034 2035 if (present(ffnl)) then 2036 if (allocated(bandfft_kpt_out%ffnl_gather)) then 2037 ABI_FREE(bandfft_kpt_out%ffnl_gather) 2038 end if 2039 if (allocated(ffnl)) then 2040 is1=size(ffnl,1) 2041 is2=size(ffnl,2) 2042 is3=size(ffnl,3) 2043 is4=size(ffnl,4) 2044 ABI_MALLOC(bandfft_kpt_out%ffnl_gather,(is1,is2,is3,is4)) 2045 bandfft_kpt_out%ffnl_gather(:,:,:,:)=ffnl(:,:,:,:) 2046 ABI_FREE(ffnl) 2047 bandfft_kpt_out%flag2_is_allocated=1 2048 end if 2049 end if 2050 if (present(ph3d)) then 2051 if (allocated(bandfft_kpt_out%ph3d_gather)) then 2052 ABI_FREE(bandfft_kpt_out%ph3d_gather) 2053 end if 2054 if (allocated(ph3d)) then 2055 is1=size(ph3d,1) 2056 is2=size(ph3d,2) 2057 is3=size(ph3d,3) 2058 ABI_MALLOC(bandfft_kpt_out%ph3d_gather,(is1,is2,is3)) 2059 bandfft_kpt_out%ph3d_gather(:,:,:)=ph3d(:,:,:) 2060 ABI_FREE(ph3d) 2061 bandfft_kpt_out%flag2_is_allocated=1 2062 end if 2063 end if 2064 if (present(kpg)) then 2065 if (allocated(bandfft_kpt_out%kpg_k_gather)) then 2066 ABI_FREE(bandfft_kpt_out%kpg_k_gather) 2067 end if 2068 if (allocated(kpg)) then 2069 is1=size(kpg,1) 2070 is2=size(kpg,2) 2071 ABI_MALLOC(bandfft_kpt_out%kpg_k_gather,(is1,is2)) 2072 bandfft_kpt_out%kpg_k_gather(:,:)=kpg(:,:) 2073 ABI_FREE(kpg) 2074 bandfft_kpt_out%flag2_is_allocated=1 2075 end if 2076 end if 2077 if (present(kinpw)) then 2078 2079 #if defined HAVE_GPU && defined HAVE_YAKL 2080 if (associated(bandfft_kpt_out%kinpw_gather)) then 2081 if(bandfft_kpt_out%gpu_option==ABI_GPU_KOKKOS) then 2082 ABI_FREE_MANAGED(bandfft_kpt_out%kinpw_gather) 2083 else 2084 ABI_FREE(bandfft_kpt_out%kinpw_gather) 2085 end if 2086 end if 2087 if (allocated(kinpw)) then 2088 is1=size(kinpw,1) 2089 if(bandfft_kpt_out%gpu_option==ABI_GPU_KOKKOS) then 2090 ABI_MALLOC_MANAGED(bandfft_kpt_out%kinpw_gather,(/is1/)) 2091 else 2092 ABI_MALLOC(bandfft_kpt_out%kinpw_gather,(is1)) 2093 end if 2094 bandfft_kpt_out%kinpw_gather(:)=kinpw(:) 2095 ABI_FREE(kinpw) 2096 end if 2097 #else 2098 if (allocated(bandfft_kpt_out%kinpw_gather)) then 2099 ABI_FREE(bandfft_kpt_out%kinpw_gather) 2100 end if 2101 if (allocated(kinpw)) then 2102 is1=size(kinpw,1) 2103 ABI_MALLOC(bandfft_kpt_out%kinpw_gather,(is1)) 2104 bandfft_kpt_out%kinpw_gather(:)=kinpw(:) 2105 ABI_FREE(kinpw) 2106 end if 2107 #endif 2108 2109 end if 2110 2111 end subroutine bandfft_kpt_restoretabs
m_bandfft_kpt/bandfft_kpt_savetabs [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_savetabs
FUNCTION
Save some (kpt-dependent) tabs from a bandfft_kpt datastructure
INPUTS
bandfft_kpt_in=<type(bandfft_kpt)>=bandfft_kpt datastructure
OUTPUT
ffnl(:,:,:,:)=nonlocal form factors on basis sphere ph3d(:,:,:)=3-dim structure factors, for each atom and plane wave kpg(:,:)=k+G vector for a given k point kinpw(:)=kinetic energy for each plane wave (Hartree)
SOURCE
1932 subroutine bandfft_kpt_savetabs(bandfft_kpt_in,ffnl,ph3d,kpg,kinpw) 1933 1934 !Arguments ------------------------------- 1935 type(bandfft_kpt_type),intent(inout) :: bandfft_kpt_in 1936 real(dp),intent(inout),allocatable,optional :: ffnl(:,:,:,:),ph3d(:,:,:),kpg(:,:),kinpw(:) 1937 !Local variables------------------------------- 1938 integer :: is1,is2,is3,is4 1939 1940 ! ********************************************************************* 1941 1942 if (present(ffnl)) then 1943 if (allocated(ffnl)) then 1944 ABI_FREE(ffnl) 1945 end if 1946 if (allocated(bandfft_kpt_in%ffnl_gather)) then 1947 is1=size(bandfft_kpt_in%ffnl_gather,1) 1948 is2=size(bandfft_kpt_in%ffnl_gather,2) 1949 is3=size(bandfft_kpt_in%ffnl_gather,3) 1950 is4=size(bandfft_kpt_in%ffnl_gather,4) 1951 ABI_MALLOC(ffnl,(is1,is2,is3,is4)) 1952 ffnl(:,:,:,:)=bandfft_kpt_in%ffnl_gather(:,:,:,:) 1953 end if 1954 end if 1955 if (present(ph3d)) then 1956 if (allocated(ph3d)) then 1957 ABI_FREE(ph3d) 1958 end if 1959 if (allocated(bandfft_kpt_in%ph3d_gather)) then 1960 is1=size(bandfft_kpt_in%ph3d_gather,1) 1961 is2=size(bandfft_kpt_in%ph3d_gather,2) 1962 is3=size(bandfft_kpt_in%ph3d_gather,3) 1963 ABI_MALLOC(ph3d,(is1,is2,is3)) 1964 ph3d(:,:,:)=bandfft_kpt_in%ph3d_gather(:,:,:) 1965 end if 1966 end if 1967 if (present(kpg)) then 1968 if (allocated(kpg)) then 1969 ABI_FREE(kpg) 1970 end if 1971 if (allocated(bandfft_kpt_in%kpg_k_gather)) then 1972 is1=size(bandfft_kpt_in%kpg_k_gather,1) 1973 is2=size(bandfft_kpt_in%kpg_k_gather,2) 1974 ABI_MALLOC(kpg,(is1,is2)) 1975 kpg(:,:)=bandfft_kpt_in%kpg_k_gather(:,:) 1976 end if 1977 end if 1978 if (present(kinpw)) then 1979 if (allocated(kinpw)) then 1980 ABI_FREE(kinpw) 1981 end if 1982 1983 #if defined HAVE_GPU && defined HAVE_YAKL 1984 if (associated(bandfft_kpt_in%kinpw_gather)) then 1985 is1=size(bandfft_kpt_in%kinpw_gather,1) 1986 ABI_MALLOC(kinpw,(is1)) 1987 kinpw(:)=bandfft_kpt_in%kinpw_gather(:) 1988 end if 1989 #else 1990 if (allocated(bandfft_kpt_in%kinpw_gather)) then 1991 is1=size(bandfft_kpt_in%kinpw_gather,1) 1992 ABI_MALLOC(kinpw,(is1)) 1993 kinpw(:)=bandfft_kpt_in%kinpw_gather(:) 1994 end if 1995 #endif 1996 1997 end if 1998 1999 end subroutine bandfft_kpt_savetabs
m_bandfft_kpt/bandfft_kpt_set_ikpt [ Functions ]
[ Top ] [ m_bandfft_kpt ] [ Functions ]
NAME
bandfft_kpt_set_ikpt
FUNCTION
Set the value of bandfft_kpt_current_ikpt variable, i.e. current index of bandfft_kpt loaded in memory INPUT ikpt=index of k-point (in the global array dtset%kpt) mpi_enreg= information about MPI parallelization
OUTPUT
bandfft_kpt_current_ikpt value changed
SOURCE
2133 subroutine bandfft_kpt_set_ikpt(ikpt,mpi_enreg) 2134 2135 !Arguments ------------------------------- 2136 integer,intent(in) :: ikpt 2137 type(MPI_type),intent(inout) :: mpi_enreg 2138 !Local variables------------------------------- 2139 2140 ! ********************************************************************* 2141 2142 if (ikpt>0) then 2143 bandfft_kpt_current_ikpt=mpi_enreg%my_kpttab(ikpt) 2144 else 2145 bandfft_kpt_current_ikpt=-1 2146 end if 2147 2148 end subroutine bandfft_kpt_set_ikpt
m_bandfft_kpt/bandfft_kpt_type [ Types ]
[ Top ] [ m_bandfft_kpt ] [ Types ]
NAME
bandfft_kpt_type
FUNCTION
The bandfft_kpt_type structured datatype gather different information about the triple band-fft-kpt parallelisation: tabs which are distributed over all the three dimensions and stored during the calculation, dimensions of messages exchange during the calculations... i.e.: all the information which were spread over the entire code before and recomputed at each iline, istep or itime STEP with a large probability to make a mistake.
SOURCE
76 type, public :: bandfft_kpt_type 77 78 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 79 ! declared in another part of ABINIT, that might need to take into account your modification. 80 81 integer :: flag1_is_allocated ! determine if the following data are allocated or not 82 integer :: npw_tot ! array holding the total number of plane waves for each k point 83 integer :: ndatarecv ! total number of values received by the processor and sent 84 ! by the other processors band 85 86 #if defined HAVE_GPU && defined HAVE_YAKL 87 integer(c_int32_t), ABI_CONTIGUOUS pointer :: kg_k_gather(:,:) => null() 88 #else 89 integer, allocatable :: kg_k_gather(:,:) ! planewave coordinates 90 ! (of the processor + sent by other processors band) 91 #endif 92 93 integer, allocatable :: recvcounts(:) ! number of values received by the processor from each processor band 94 integer, allocatable :: sendcounts(:) ! number of values sent by the processor to each processor band 95 integer, allocatable :: rdispls (:) ! positions of values received by the processor from each processor band 96 integer, allocatable :: sdispls (:) ! postions of values sent by the processor to each processor band 97 integer, allocatable :: gbound(:,:) ! sphere boundary info: gbound(2*mgfft+8,2) 98 99 integer :: flag2_is_allocated ! determine if the following data are allocated or not 100 real(dp), allocatable :: ffnl_gather(:,:,:,:) ! ffnl tab (of the processor + sent by other processors band) 101 102 #if defined HAVE_GPU && defined HAVE_YAKL 103 real(c_double), ABI_CONTIGUOUS pointer :: kinpw_gather(:) => null() ! kinpw tab (of the processor + sent by other processors band) 104 #else 105 real(dp), allocatable :: kinpw_gather(:) ! kinpw tab (of the processor + sent by other processors band) 106 #endif 107 108 real(dp), allocatable :: ph3d_gather(:,:,:) ! ph3d tab (of the processor + sent by other processors band) 109 real(dp), allocatable :: kpg_k_gather(:,:) ! kpg_k tab (of the processor + sent by other processors band) 110 111 112 integer :: flag3_is_allocated ! determine if the following data are allocated or not 113 integer :: istwf_k ! input option parameter that describes the storage of wfs 114 integer :: idatarecv0 ! position of the planewave coordinates (0,0,0) 115 integer :: ndatarecv_tot ! total number of received values by the processor 116 ! (ndatarecv + number of received opposited planewave coordinates) 117 integer :: ndatasend_sym ! number of sent values to the processors fft to create opposited 118 ! planewave coordinates 119 integer, allocatable :: kg_k_gather_sym(:,:) ! planewave coordinates 120 ! (kg_k_gather + opposited planewave coordinates sent by the processors fft) 121 integer, allocatable :: rdispls_sym(:) ! positions of values received by the processor from each processor fft 122 integer, allocatable :: recvcounts_sym(:) ! number of values received by the processor from each processor fft 123 integer, allocatable :: recvcounts_sym_tot(:) ! number of values received by each processor from the other processors fft 124 integer, allocatable :: sdispls_sym(:) ! postions of values sent by the processor to each processor fft 125 integer, allocatable :: sendcounts_sym(:) ! number of values sent by the processor to each processor fft 126 integer, allocatable :: sendcounts_sym_all(:) ! number of values sent by each processor to the other processors fft 127 integer, allocatable :: tab_proc(:) ! positions of opposited planewave coordinates in the list of the processors fft 128 129 logical :: have_to_reequilibrate ! indicates weather we will have to reequilibrate and allocate all these stuff 130 integer :: npw_fft ! Number of plane waves during fft step 131 integer, allocatable :: indices_pw_fft(:) ! Indices for sorting pw like pw 132 integer, allocatable :: sendcount_fft (:) ! Number of pw to send to others proc fft 133 integer, allocatable :: senddisp_fft(:) ! Positions for sending 134 integer, allocatable :: recvcount_fft(:) ! Number of pw to receive from others proc fft 135 integer, allocatable :: recvdisp_fft(:) ! Positions for receiving 136 integer, allocatable :: kg_k_fft(:,:) ! planewaves coordinates 137 integer :: gpu_option ! if this structure will be used with GPU 138 139 end type bandfft_kpt_type