Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : module m_wann_write_amn
8 : use m_juDFT
9 : contains
10 4 : subroutine wann_write_amn(
11 : > fmpi_comm,l_p0,filename,title,
12 : > nbnd,fullnkpts,nwfs,
13 : > irank,isize,l_freeformat,l_unfname,
14 4 : < amn,l_unformatted)
15 : c**********************************************************
16 : c This subroutine is used to write several matrices to
17 : c files: WF1.mmn, WF1.amn, etc. The corresponding
18 : c filename has to be provided as input. To be concrete
19 : c all explanations given in the following refer to
20 : c WF1.amn/WF2.amn.
21 : c
22 : c MPI-Version: Collect the contributions to the matrix
23 : c A^{k}_{mn} from the various processors.
24 : c
25 : c Write the matrix A^{k}_{mn} to file WF1.amn/WF2.amn
26 : c
27 : c Frank Freimuth
28 : c**********************************************************
29 : #ifdef CPP_MPI
30 : use mpi
31 : #endif
32 :
33 : implicit none
34 : integer, intent(in) :: fmpi_comm
35 : logical, intent(in) :: l_p0,l_unformatted
36 : character, intent(in) :: filename*(*)
37 : character, intent(in) :: title*(*)
38 :
39 : integer, intent(in) :: nbnd
40 : integer, intent(in) :: fullnkpts
41 : integer, intent(in) :: nwfs
42 :
43 : integer, intent(in) :: irank,isize
44 : logical, intent(in) :: l_freeformat
45 : logical, intent(in) :: l_unfname
46 :
47 : complex, intent(inout) :: amn(:,:,:)
48 :
49 : integer :: ikpt,nwf,i
50 : integer :: cpu_index
51 : #ifdef CPP_MPI
52 : integer :: ierr(3)
53 : integer :: stt(MPI_STATUS_SIZE)
54 : #endif
55 :
56 4 : call timestart("wann_write_amn")
57 : #ifdef CPP_MPI
58 : c******************************************************
59 : c Collect contributions to the amn matrix from the
60 : c various processors.
61 : c******************************************************
62 :
63 4 : if(isize.ne.1)then
64 36 : do ikpt=1,fullnkpts
65 32 : if(l_p0)then
66 32 : do cpu_index=1,isize-1
67 32 : if(mod(ikpt-1,isize).eq.cpu_index)then
68 : call MPI_RECV(
69 : & amn(1:nbnd,1:nwfs,ikpt),nbnd*nwfs,
70 : & MPI_DOUBLE_COMPLEX,cpu_index,
71 8 : & ikpt,fmpi_comm,stt,ierr(1))
72 : endif !processors
73 : enddo !cpu_index
74 : else
75 16 : if(mod(ikpt-1,isize).eq.irank)then
76 : call MPI_SEND(
77 : & amn(1:nbnd,1:nwfs,ikpt),nbnd*nwfs,
78 : & MPI_DOUBLE_COMPLEX,0,
79 8 : & ikpt,fmpi_comm,ierr(1))
80 : endif !processors
81 : endif ! l_p0
82 36 : call MPI_BARRIER(fmpi_comm,ierr(1))
83 : enddo !ikpt
84 : endif !isize
85 : #endif
86 :
87 4 : if(l_p0)then
88 2 : if(l_unformatted)then
89 0 : if(l_unfname)then
90 0 : open (305,file=trim(filename)//'_unf',form='unformatted')
91 : else
92 0 : open (305,file=filename,form='unformatted')
93 : endif
94 2 : elseif(l_freeformat)then
95 0 : open(305,file=filename,recl=1000)
96 : else
97 2 : open (305,file=filename)
98 : endif
99 :
100 2 : if(l_unformatted)then
101 0 : write(305)nbnd,fullnkpts,nwfs
102 0 : write(305)amn
103 : else
104 2 : write (305,*)title
105 2 : write (305,'(i5,i7,i5)') nbnd,fullnkpts,nwfs
106 2 : if(l_freeformat)then
107 0 : do ikpt = 1,fullnkpts
108 0 : do nwf = 1,nwfs
109 0 : do i = 1,nbnd
110 0 : write (305,*) i,nwf,ikpt,
111 0 : & real(amn(i,nwf,ikpt)),aimag(amn(i,nwf,ikpt))
112 : enddo !i
113 : enddo !nwf
114 : enddo !ikpt
115 : else
116 18 : do ikpt = 1,fullnkpts
117 146 : do nwf = 1,nwfs
118 1168 : do i = 1,nbnd
119 1024 : write (305,'(i5,i5,i7,3x,2f18.12)') i,nwf,ikpt,
120 2176 : & real(amn(i,nwf,ikpt)),aimag(amn(i,nwf,ikpt))
121 : enddo !i
122 : enddo !nwf
123 : enddo !ikpt
124 : endif
125 : endif
126 2 : close(305)
127 : endif
128 4 : call timestop("wann_write_amn")
129 :
130 4 : end subroutine wann_write_amn
131 : end module m_wann_write_amn
|