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_matrix6
8 : use m_juDFT
9 : contains
10 0 : subroutine wann_write_matrix6(
11 : > fmpi_comm,l_p0,filename,title,
12 : > spin1,spin2,num_bands1,num_bands2,
13 : > num_dims,fullnkpts,
14 : > irank,isize,l_unformatted,
15 0 : < matrix6)
16 : c*************************************************************
17 : c This subroutine is used to write several matrices to
18 : c files. The corresponding
19 : c filename has to be provided as input.
20 : c
21 : c MPI-Version: Collect the contributions to the matrix
22 : c from the various processors.
23 : c
24 : c Frank Freimuth
25 : c*************************************************************
26 :
27 : USE m_constants
28 : #ifdef CPP_MPI
29 : USE mpi
30 : #endif
31 :
32 : implicit none
33 :
34 : integer, intent(in) :: fmpi_comm
35 : logical, intent(in) :: l_p0
36 : character, intent(in) :: filename*(*)
37 : character, intent(in) :: title*(*)
38 :
39 : integer, intent(in) :: spin1
40 : integer, intent(in) :: spin2
41 : integer, intent(in) :: num_bands1
42 : integer, intent(in) :: num_bands2
43 :
44 : integer, intent(in) :: num_dims
45 : integer, intent(in) :: fullnkpts
46 :
47 : integer, intent(in) :: irank,isize
48 : logical, intent(in) :: l_unformatted
49 :
50 : complex, intent(inout) :: matrix6(:,:,:,:,:,:)
51 :
52 : integer :: ikpt,i,j,ii,jj
53 : integer :: cpu_index,dir
54 : #ifdef CPP_MPI
55 : integer :: ierr(3)
56 : integer :: stt(MPI_STATUS_SIZE)
57 : #endif
58 :
59 : #ifdef CPP_MPI
60 : c**********************************************************
61 : c Collect contributions to the matrix6 matrix from the
62 : c various processors.
63 : c**********************************************************
64 :
65 0 : call timestart("wann_write_matrix6")
66 :
67 0 : if(isize.ne.1)then
68 0 : do ikpt=1,fullnkpts
69 0 : if(l_p0)then
70 0 : do cpu_index=1,isize-1
71 0 : if(mod(ikpt-1,isize).eq.cpu_index)then
72 : call MPI_RECV(
73 : & matrix6(1:spin1,1:spin2,1:num_bands1,
74 : & 1:num_bands2,1:num_dims,ikpt),
75 : & spin1*spin2*num_bands1*num_bands2*num_dims,
76 : & MPI_DOUBLE_COMPLEX,cpu_index,
77 0 : & ikpt,fmpi_comm,stt,ierr(1))
78 : endif !processors
79 : enddo !cpu_index
80 : else
81 0 : if(mod(ikpt-1,isize).eq.irank)then
82 : call MPI_SEND(
83 : & matrix6(1:spin1,1:spin2,1:num_bands1,
84 : & 1:num_bands2,1:num_dims,ikpt),
85 : & spin1*spin2*num_bands1*num_bands2*num_dims,
86 : & MPI_DOUBLE_COMPLEX,0,
87 0 : & ikpt,fmpi_comm,ierr(1))
88 : endif !processors
89 : endif ! l_p0
90 0 : call MPI_BARRIER(fmpi_comm,ierr(1))
91 : enddo !ikpt
92 : endif !isize
93 : #endif
94 :
95 0 : write(oUnit,*)"wann_write_matrix6"
96 :
97 0 : if(l_p0)then
98 0 : if(l_unformatted)then
99 0 : open(305,file=trim(filename)//'_unf',form='unformatted')
100 0 : write(305)spin2,spin1,num_bands1,num_bands2,
101 0 : & num_dims,fullnkpts
102 0 : do ikpt=1,fullnkpts
103 0 : write(305)matrix6(:,:,:,:,:,ikpt)
104 : enddo
105 : else
106 0 : open (305,file=filename)
107 0 : write (305,*)title
108 0 : write (305,'(6i5)') spin2,spin1,num_bands1,num_bands2,
109 0 : & num_dims,fullnkpts
110 0 : do ikpt=1,fullnkpts
111 0 : do dir=1,num_dims
112 0 : do i = 1,num_bands2
113 0 : do j = 1,num_bands1
114 0 : do ii=1,spin1
115 0 : do jj=1,spin2
116 0 : write (305,'(6i5,3x,2f18.12)') jj,ii,j,i,dir,ikpt,
117 0 : & real(matrix6(jj,ii,j,i,dir,ikpt)),
118 0 : & aimag(matrix6(jj,ii,j,i,dir,ikpt))
119 : enddo !jj
120 : enddo !ii
121 : enddo !j
122 : enddo !i
123 : enddo !dir
124 : enddo !ikpt
125 : endif
126 0 : close(305)
127 : endif
128 0 : call timestop("wann_write_matrix6")
129 0 : end subroutine wann_write_matrix6
130 : end module m_wann_write_matrix6
|