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_mmnk
8 : #ifdef CPP_MPI
9 : use mpi
10 : #endif
11 : contains
12 2 : subroutine wann_write_mmnk(
13 : > mpi_communicatior,jspin2,l_p0,fullnkpts,nntot,
14 4 : > wann,maptopair,pair_to_do,nbnd,bpt,gb,
15 : > isize,irank,fending,
16 2 : < mmnk,l_unformatted)
17 : c**********************************************************
18 : c MPI-Version: Collect the contributions to the matrix
19 : c M^{k,b}_{mn} from the various processors.
20 : c
21 : c Symmetry: Compose the M^{k,b}_{mn} matrix for the
22 : c full BZ from the pieces of the IBZ.
23 : c
24 : c Write the matrix M^{k,b}_{mn} to file WF1.mmn/WF2.mmn
25 : c Frank Freimuth
26 : c**********************************************************
27 : use m_types
28 : use m_juDFT
29 : implicit none
30 : integer, intent(in) :: jspin2,mpi_communicatior
31 : logical, intent(in) :: l_p0,l_unformatted
32 : integer, intent(in) :: fullnkpts
33 : integer, intent(in) :: nntot
34 : type(t_wann),intent(in) :: wann
35 :
36 : integer, intent(in) :: maptopair(:,:,:) !maptopair(3,fullnkpts,nntot)
37 : integer, intent(in) :: pair_to_do(:,:) !pair_to_do(fullnkpts,nntot)
38 : integer, intent(in) :: nbnd
39 : integer, intent(in) :: bpt(:,:)
40 : integer, intent(in) :: gb(:,:,:)
41 :
42 : integer, intent(in) :: isize,irank
43 :
44 : CHARACTER(len=12), INTENT(IN) :: fending !for file ending !QPOINTS
45 :
46 : complex, intent(inout) :: mmnk(:,:,:,:) !mmnk(nbnd,nbnd,nntot,fullnkpts)
47 :
48 : integer :: ikpt,i,j
49 : integer :: ikpt_b
50 : character(len=3) :: spin12(2)
51 : integer :: cpu_index
52 : data spin12/'WF1' , 'WF2'/
53 :
54 : #ifdef CPP_MPI
55 : integer :: ierr
56 : integer :: stt(MPI_STATUS_SIZE)
57 : #endif
58 :
59 2 : call timestart("wann_write_mmnk")
60 :
61 : #ifdef CPP_MPI
62 : c******************************************************
63 : c Collect contributions to the mmnk matrix from the
64 : c various processors.
65 : c******************************************************
66 2 : if(isize.ne.1)then
67 18 : do ikpt=1,fullnkpts
68 16 : if(l_p0)then
69 16 : do cpu_index=1,isize-1
70 16 : if(mod(ikpt-1,isize).eq.cpu_index)then
71 36 : do ikpt_b=1,nntot !nearest neighbors
72 36 : if(pair_to_do(ikpt,ikpt_b).ne.0)then
73 : call MPI_RECV(
74 : & mmnk(1:nbnd,1:nbnd,ikpt_b,ikpt),nbnd*nbnd,
75 : & MPI_DOUBLE_COMPLEX,cpu_index,5*fullnkpts+
76 12 : & pair_to_do(ikpt,ikpt_b),mpi_communicatior,stt,ierr)
77 :
78 : endif !pairtodo
79 : enddo !nearest neighbors
80 : endif !processors
81 : enddo !cpu_index
82 : else
83 8 : if(mod(ikpt-1,isize).eq.irank)then
84 36 : do ikpt_b=1,nntot !loop over nearest neighbors
85 36 : if(pair_to_do(ikpt,ikpt_b).ne.0)then
86 : call MPI_SEND(
87 : & mmnk(1:nbnd,1:nbnd,ikpt_b,ikpt),
88 : & nbnd*nbnd,MPI_DOUBLE_COMPLEX,0,5*fullnkpts+
89 12 : & pair_to_do(ikpt,ikpt_b),mpi_communicatior,ierr)
90 : endif !pairtodo
91 : enddo !loop over nearest neighbors
92 : endif !processors
93 : endif ! l_p0
94 18 : call MPI_BARRIER(mpi_communicatior,ierr)
95 : enddo !ikpt
96 : endif !isize
97 : #endif
98 :
99 : c****************************************************
100 : c Symmetry: complete the mmnk matrix.
101 : c****************************************************
102 2 : if(l_p0)then
103 9 : do ikpt=1,fullnkpts
104 73 : do ikpt_b=1,nntot
105 72 : if(pair_to_do(ikpt,ikpt_b).eq.0)then
106 32 : if(maptopair(3,ikpt,ikpt_b).eq.1)then !conjugation selection
107 : mmnk(:,:,ikpt_b,ikpt)=conjg(transpose(mmnk(:,:,
108 4640 : & maptopair(2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b))))
109 0 : elseif(maptopair(3,ikpt,ikpt_b).eq.2)then !rotation
110 : mmnk(:,:,ikpt_b,ikpt)=mmnk(:,:,maptopair
111 0 : & (2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b))
112 0 : elseif(maptopair(3,ikpt,ikpt_b).eq.3)then !rotation&reflection
113 : mmnk(:,:,ikpt_b,ikpt)=conjg( mmnk(:,:,maptopair
114 0 : & (2,ikpt,ikpt_b),maptopair(1,ikpt,ikpt_b)) )
115 : else !something wrong
116 0 : call juDFT_error('maptopair')
117 : endif!maptopair
118 : endif!pairtodo
119 : enddo!ikpt_b
120 : enddo!ikpt
121 : endif
122 :
123 : c******************************************************
124 : c Write mmnk matrix to file.
125 : c******************************************************
126 2 : if (l_p0) then
127 1 : if(.not.l_unformatted) then
128 1 : open (305,file=spin12(jspin2)//trim(fending)//'.mmn')
129 1 : write (305,*) 'Overlaps of the wavefunct. the k- and b-points'
130 1 : write (305,'(3i5)') nbnd,fullnkpts,nntot
131 9 : do ikpt = 1,fullnkpts
132 73 : do ikpt_b = 1,nntot
133 64 : write (305,'(2i5,3x,3i4)') ikpt,bpt(ikpt_b,ikpt),
134 128 : & gb(1:3,ikpt_b,ikpt)
135 584 : do i = 1,nbnd
136 4672 : do j = 1,nbnd
137 : c write (305,'(2f18.12)')
138 : write (305,'(2f24.18)')
139 4608 : & real(mmnk(j,i,ikpt_b,ikpt)),-aimag(mmnk(j,i,ikpt_b,ikpt))
140 : enddo
141 : enddo
142 : enddo
143 : enddo !ikpt
144 1 : close (305)
145 : else
146 : open (305,file=spin12(jspin2)//trim(fending)//'.mmn',
147 0 : > form='unformatted')
148 0 : write (305) nbnd,fullnkpts,nntot
149 0 : write (305) bpt,gb
150 0 : write (305) conjg(mmnk)
151 0 : close (305)
152 : endif
153 : endif !l_p0
154 2 : call timestop("wann_write_mmnk")
155 2 : end subroutine wann_write_mmnk
156 : end module m_wann_write_mmnk
|