wann_write_amn.F 3.81 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_amn
      contains
      subroutine wann_write_amn(
     >               mpi_comm,l_p0,filename,title,
     >               nbnd,fullnkpts,nwfs,
Daniel Wortmann's avatar
Daniel Wortmann committed
12 13
     >               irank,isize,l_freeformat,
     <               amn,l_unformatted)
14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29
c**********************************************************
c     This subroutine is used to write several matrices to 
c     files: WF1.mmn, WF1.amn, etc. The corresponding 
c     filename has to be provided as input. To be concrete
c     all explanations given in the following refer to
c     WF1.amn/WF2.amn.
c
c     MPI-Version: Collect the contributions to the matrix
c     A^{k}_{mn} from the various processors.
c
c     Write the matrix A^{k}_{mn} to file WF1.amn/WF2.amn
c
c     Frank Freimuth
c**********************************************************
      implicit none
      integer, intent(in)    :: mpi_comm
Daniel Wortmann's avatar
Daniel Wortmann committed
30
      logical, intent(in)    :: l_p0,l_unformatted
31 32 33 34 35 36 37 38 39 40 41 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
      character, intent(in)  :: filename*(*)
      character, intent(in)  :: title*(*)

      integer, intent(in)    :: nbnd
      integer, intent(in)    :: fullnkpts
      integer, intent(in)    :: nwfs

      integer, intent(in)    :: irank,isize
      logical, intent(in)    :: l_freeformat

      complex, intent(inout) :: amn(:,:,:)

      integer :: ikpt,nwf,i
      integer :: cpu_index
#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 amn 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
           call MPI_RECV(
     &             amn(1:nbnd,1:nwfs,ikpt),nbnd*nwfs,
     &             CPP_MPI_COMPLEX,cpu_index,
     &             ikpt,mpi_comm,stt,ierr)
          endif !processors
         enddo !cpu_index
        else
         if(mod(ikpt-1,isize).eq.irank)then
           call MPI_SEND(
     &             amn(1:nbnd,1:nwfs,ikpt),nbnd*nwfs,
     &             CPP_MPI_COMPLEX,0,
     &             ikpt,mpi_comm,ierr)
         endif !processors
        endif ! l_p0
        call MPI_BARRIER(mpi_comm,ierr)
       enddo !ikpt 
      endif !isize
#endif
Daniel Wortmann's avatar
Daniel Wortmann committed
80

81 82
      if(l_p0)then
       if(l_unformatted)then
Daniel Wortmann's avatar
Daniel Wortmann committed
83 84
          open (305,file=filename,form='unformatted')
       elseif(l_freeformat)then
85
          open(305,file=filename,recl=1000)
Daniel Wortmann's avatar
Daniel Wortmann committed
86
       else
87
          open (305,file=filename)
Daniel Wortmann's avatar
Daniel Wortmann committed
88 89 90 91 92 93
       endif
       
       if(l_unformatted)then
        write(305)nbnd,fullnkpts,nwfs
        write(305)amn
       else
94
        write (305,*)title
Daniel Wortmann's avatar
Daniel Wortmann committed
95
        write (305,'(i5,i7,i5)') nbnd,fullnkpts,nwfs
96 97 98 99 100 101 102 103 104 105 106 107 108
        if(l_freeformat)then
         do ikpt = 1,fullnkpts
          do nwf = 1,nwfs
           do i = 1,nbnd
            write (305,*) i,nwf,ikpt,
     &              real(amn(i,nwf,ikpt)),aimag(amn(i,nwf,ikpt))
           enddo !i
          enddo !nwf
         enddo !ikpt
        else
         do ikpt = 1,fullnkpts
          do nwf = 1,nwfs
           do i = 1,nbnd
Daniel Wortmann's avatar
Daniel Wortmann committed
109
            write (305,'(i5,i5,i7,3x,2f18.12)') i,nwf,ikpt,
110 111 112 113
     &              real(amn(i,nwf,ikpt)),aimag(amn(i,nwf,ikpt))
           enddo !i
          enddo !nwf
         enddo !ikpt
Daniel Wortmann's avatar
Daniel Wortmann committed
114 115
        endif 
       endif
116 117 118 119 120
       close(305)
      endif

      end subroutine wann_write_amn
      end module m_wann_write_amn