wann_write_mmnk.F 5.55 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7 8 9 10 11
      module m_wann_write_mmnk
      contains
      subroutine wann_write_mmnk(
     >               mpi_comm,jspin2,l_p0,fullnkpts,nntot,wann,
     >               maptopair,pair_to_do,nbnd,bpt,gb,
Daniel Wortmann's avatar
Daniel Wortmann committed
12 13
     >               isize,irank,fending,
     <               mmnk,l_unformatted)
14 15 16 17 18 19 20 21 22 23 24 25 26 27
c**********************************************************
c     MPI-Version: Collect the contributions to the matrix
c     M^{k,b}_{mn} from the various processors.
c
c     Symmetry: Compose the M^{k,b}_{mn} matrix for the
c     full BZ from the pieces of the IBZ.
c
c     Write the matrix M^{k,b}_{mn} to file WF1.mmn/WF2.mmn
c     Frank Freimuth
c**********************************************************
      use m_types
      use m_juDFT
      implicit none
      integer, intent(in)     :: jspin2,mpi_comm
Daniel Wortmann's avatar
Daniel Wortmann committed
28
      logical, intent(in)     :: l_p0,l_unformatted
29 30 31 32 33 34 35 36 37 38 39
      integer, intent(in)     :: fullnkpts
      integer, intent(in)     :: nntot
      type(t_wann),intent(in) :: wann

      integer, intent(in)     :: maptopair(:,:,:) !maptopair(3,fullnkpts,nntot)
      integer, intent(in)     :: pair_to_do(:,:)  !pair_to_do(fullnkpts,nntot)
      integer, intent(in)     :: nbnd
      integer, intent(in)     :: bpt(:,:)
      integer, intent(in)     :: gb(:,:,:)

      integer, intent(in)     :: isize,irank
Daniel Wortmann's avatar
Daniel Wortmann committed
40 41
      
      CHARACTER(len=12), INTENT(IN) :: fending  !for file ending !QPOINTS
42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123

      complex, intent(inout)  :: mmnk(:,:,:,:) !mmnk(nbnd,nbnd,nntot,fullnkpts)

      integer          :: ikpt,i,j
      integer          :: ikpt_b
      character(len=3) :: spin12(2)
      integer          :: cpu_index
      data   spin12/'WF1' , 'WF2'/

#ifdef CPP_MPI
      include 'mpif.h'
      integer :: ierr(3)
      integer :: stt(MPI_STATUS_SIZE)
#include "cpp_double.h"
#endif

#ifdef CPP_MPI
c******************************************************
c     Collect contributions to the mmnk matrix from the
c     various processors.
c******************************************************
      if(isize.ne.1)then
      do ikpt=1,fullnkpts
       if(l_p0)then
        do cpu_index=1,isize-1
         if(mod(ikpt-1,isize).eq.cpu_index)then
          do ikpt_b=1,nntot !nearest neighbors
           if(pair_to_do(ikpt,ikpt_b).ne.0)then
             call MPI_RECV(
     &               mmnk(1:nbnd,1:nbnd,ikpt_b,ikpt),nbnd*nbnd,
     &               CPP_MPI_COMPLEX,cpu_index,5*fullnkpts+
     &               pair_to_do(ikpt,ikpt_b),mpi_comm,stt,ierr)

           endif !pairtodo
          enddo !nearest neighbors
         endif !processors
        enddo !cpu_index
       else
        if(mod(ikpt-1,isize).eq.irank)then
         do ikpt_b=1,nntot !loop over nearest neighbors
          if(pair_to_do(ikpt,ikpt_b).ne.0)then
            call MPI_SEND(
     &              mmnk(1:nbnd,1:nbnd,ikpt_b,ikpt),
     &              nbnd*nbnd,CPP_MPI_COMPLEX,0,5*fullnkpts+
     &              pair_to_do(ikpt,ikpt_b),mpi_comm,ierr)
          endif !pairtodo
         enddo !loop over nearest neighbors
        endif !processors
       endif ! l_p0
       call MPI_BARRIER(mpi_comm,ierr)
      enddo !ikpt
      endif !isize
#endif

c****************************************************
c     Symmetry: complete the mmnk matrix.
c****************************************************
      if(l_p0)then 
       do ikpt=1,fullnkpts
        do ikpt_b=1,nntot
         if(pair_to_do(ikpt,ikpt_b).eq.0)then
          if(maptopair(3,ikpt,ikpt_b).eq.1)then !conjugation selection
             mmnk(:,:,ikpt_b,ikpt)=conjg(transpose(mmnk(:,:,
     &          maptopair(2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b))))
          elseif(maptopair(3,ikpt,ikpt_b).eq.2)then !rotation
             mmnk(:,:,ikpt_b,ikpt)=mmnk(:,:,maptopair
     &         (2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b))
          elseif(maptopair(3,ikpt,ikpt_b).eq.3)then !rotation&reflection
             mmnk(:,:,ikpt_b,ikpt)=conjg( mmnk(:,:,maptopair
     &         (2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b)) )
          else !something wrong
             call juDFT_error('maptopair')
          endif!maptopair
         endif!pairtodo
        enddo!ikpt_b
       enddo!ikpt
      endif

c******************************************************
c     Write mmnk matrix to file.
c******************************************************
      if (l_p0) then
Daniel Wortmann's avatar
Daniel Wortmann committed
124 125
       if(.not.l_unformatted) then
       open (305,file=spin12(jspin2)//trim(fending)//'.mmn')
126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141
       write (305,*) 'Overlaps of the wavefunct. the k- and b-points'
       write (305,'(3i5)') nbnd,fullnkpts,nntot
       do ikpt = 1,fullnkpts
        do ikpt_b = 1,nntot
         write (305,'(2i5,3x,3i4)') ikpt,bpt(ikpt_b,ikpt),
     &                                  gb(1:3,ikpt_b,ikpt)
         do i = 1,nbnd
          do j = 1,nbnd
c           write (305,'(2f18.12)')
           write (305,'(2f24.18)')
     &         real(mmnk(j,i,ikpt_b,ikpt)),-aimag(mmnk(j,i,ikpt_b,ikpt))
          enddo
         enddo
        enddo
       enddo !ikpt
       close (305)
Daniel Wortmann's avatar
Daniel Wortmann committed
142 143 144 145 146 147 148 149
       else
       open (305,file=spin12(jspin2)//trim(fending)//'.mmn',
     >       form='unformatted')
       write (305) nbnd,fullnkpts,nntot
       write (305) bpt,gb
       write (305) conjg(mmnk)
       close (305)
       endif
150 151 152 153
      endif !l_p0

      end subroutine wann_write_mmnk
      end module m_wann_write_mmnk