Line data Source code
1 : module m_types_coul
2 : use m_types_mat
3 : use m_mtir_size
4 : use m_types_fleurinput
5 : use m_judft
6 : #ifdef CPP_MPI
7 : use mpi
8 : #endif
9 : implicit none
10 : private
11 : type t_coul
12 : REAL, ALLOCATABLE :: mt1_r(:, :, :, :)
13 : REAL, ALLOCATABLE :: mt2_r(:, :, :, :), mt3_r(:, :, :)
14 : COMPLEX, ALLOCATABLE :: mt1_c(:, :, :, :)
15 : COMPLEX, ALLOCATABLE :: mt2_c(:, :, :, :), mt3_c(:, :, :)
16 : type(t_mat) :: mtir
17 : #ifdef CPP_MPI
18 : integer :: comm = MPI_COMM_NULL ! communicator for this coulomb matrix
19 : #else
20 : integer :: comm = -1
21 : #endif
22 : logical :: l_participate = .False. ! am i somehow involved with this coulomb mtx
23 : contains
24 : procedure :: alloc => t_coul_alloc
25 : procedure :: mini_alloc => t_coul_mini_alloc
26 : procedure :: free => t_coul_free
27 : procedure :: mpi_bcast => t_coul_mpi_bc
28 : end type t_coul
29 : public t_coul
30 : contains
31 :
32 36 : subroutine t_coul_mpi_bc(coul, fi, communicator, root)
33 : use m_types_fleurinput
34 : use m_types_hybmpi
35 : use m_judft
36 :
37 : implicit none
38 : class(t_coul) :: coul
39 : type(t_fleurinput), intent(in) :: fi
40 : integer, intent(in) :: root, communicator
41 : #ifdef CPP_MPI
42 : integer :: ierr
43 :
44 36 : call timestart("Bcast coulomb_mtx")
45 36 : call timestart("Bcast small stuff")
46 36 : if (fi%sym%invs) THEN
47 120 : call MPI_Bcast(coul%mt1_r, size(coul%mt1_r), MPI_DOUBLE_PRECISION, root, communicator, ierr)
48 120 : call MPI_Bcast(coul%mt2_r, size(coul%mt2_r), MPI_DOUBLE_PRECISION, root, communicator, ierr)
49 96 : call MPI_Bcast(coul%mt3_r, size(coul%mt3_r), MPI_DOUBLE_PRECISION, root, communicator, ierr)
50 : else
51 60 : call MPI_Bcast(coul%mt1_c, size(coul%mt1_c), MPI_DOUBLE_COMPLEX, root, communicator, ierr)
52 60 : call MPI_Bcast(coul%mt2_c, size(coul%mt2_c), MPI_DOUBLE_COMPLEX, root, communicator, ierr)
53 48 : call MPI_Bcast(coul%mt3_c, size(coul%mt3_c), MPI_DOUBLE_COMPLEX, root, communicator, ierr)
54 : endif
55 36 : call timestop("Bcast small stuff")
56 :
57 36 : call coul%mtir%bcast(root, communicator)
58 36 : call timestop("Bcast coulomb_mtx")
59 : #endif
60 36 : end subroutine t_coul_mpi_bc
61 :
62 0 : subroutine t_coul_free(coul)
63 : implicit none
64 : class(t_coul), intent(inout) :: coul
65 : integer :: ierr
66 :
67 0 : if (allocated(coul%mt1_r)) deallocate (coul%mt1_r)
68 0 : if (allocated(coul%mt1_c)) deallocate (coul%mt1_c)
69 0 : if (allocated(coul%mt2_r)) deallocate (coul%mt2_r)
70 0 : if (allocated(coul%mt3_r)) deallocate (coul%mt3_r)
71 0 : if (allocated(coul%mt2_c)) deallocate (coul%mt2_c)
72 0 : if (allocated(coul%mt3_c)) deallocate (coul%mt3_c)
73 0 : call coul%mtir%free()
74 :
75 : #ifdef CPP_MPI
76 0 : if (coul%comm /= MPI_COMM_NULL) call MPI_comm_free(coul%comm, ierr)
77 : #endif
78 0 : end subroutine t_coul_free
79 :
80 36 : subroutine t_coul_alloc(coul, fi, num_radbasfn, n_g, ikpt, l_print)
81 : implicit NONE
82 : class(t_coul), intent(inout) :: coul
83 : type(t_fleurinput), intent(in) :: fi
84 : integer, intent(in) :: num_radbasfn(:, :), n_g(:), ikpt
85 : logical, intent(in), optional :: l_print
86 : integer :: info, isize
87 :
88 36 : isize = mtir_size(fi, n_g, ikpt)
89 :
90 36 : if (present(l_print)) then
91 36 : if (l_print) then
92 0 : write (*, *) "Coulomb dimensions:"
93 0 : write (*, *) "real:", fi%sym%invs
94 : write (*, *) "mt1 -> ["//int2str(maxval(num_radbasfn) - 1)// &
95 : ", "//int2str(maxval(num_radbasfn) - 1)// &
96 : ", "//int2str(maxval(fi%hybinp%lcutm1) + 1)// &
97 0 : ", "//int2str(fi%atoms%ntype)//"]"
98 : write (*, *) "mt2 -> ["//int2str(maxval(num_radbasfn) - 1)// &
99 : ", "//int2str(2*maxval(fi%hybinp%lcutm1) + 1)// &
100 : ", "//int2str(maxval(fi%hybinp%lcutm1) + 2)// &
101 0 : ", "//int2str(fi%atoms%nat)//"]"
102 : write (*, *) "mt3 -> ["//int2str(maxval(num_radbasfn) - 1)// &
103 : ", "//int2str(fi%atoms%nat)// &
104 0 : ", "//int2str(fi%atoms%nat)//"]"
105 : write (*, *) "mtir-> ["//int2str(isize)// &
106 0 : ", "//int2str(isize)//"]"
107 : endif
108 : endif
109 :
110 36 : if (fi%sym%invs) THEN
111 24 : if (.not. allocated(coul%mt1_r)) then
112 : allocate (coul%mt1_r(maxval(num_radbasfn) - 1, &
113 : maxval(num_radbasfn) - 1, &
114 198 : 0:maxval(fi%hybinp%lcutm1), fi%atoms%ntype), stat=info)
115 12 : if (info /= 0) call judft_error("Can't allocate coul%mt1_r")
116 : endif
117 24 : if (.not. allocated(coul%mt2_r)) then
118 : allocate (coul%mt2_r(maxval(num_radbasfn) - 1, &
119 : -maxval(fi%hybinp%lcutm1):maxval(fi%hybinp%lcutm1), &
120 216 : 0:maxval(fi%hybinp%lcutm1) + 1, fi%atoms%nat), stat=info)
121 12 : if (info /= 0) call judft_error("Can't allocate coul%mt2_r")
122 : endif
123 :
124 24 : if (.not. allocated(coul%mt3_r)) then
125 168 : allocate (coul%mt3_r(maxval(num_radbasfn) - 1, fi%atoms%nat, fi%atoms%nat), stat=info)
126 12 : if (info /= 0) call judft_error("Can't allocate coul%mt3_r")
127 : endif
128 24 : if (.not. coul%mtir%allocated()) call coul%mtir%alloc(.True., isize, isize)
129 : else
130 12 : if (.not. allocated(coul%mt1_c)) then
131 : allocate (coul%mt1_c(maxval(num_radbasfn) - 1, &
132 : maxval(num_radbasfn) - 1, &
133 120 : 0:maxval(fi%hybinp%lcutm1), fi%atoms%ntype), stat=info)
134 6 : if (info /= 0) call judft_error("Can't allocate coul%mt1_c")
135 : endif
136 12 : if (.not. allocated(coul%mt2_c)) then
137 : allocate (coul%mt2_c(maxval(num_radbasfn) - 1, &
138 : -maxval(fi%hybinp%lcutm1):maxval(fi%hybinp%lcutm1), &
139 132 : 0:maxval(fi%hybinp%lcutm1) + 1, fi%atoms%nat), stat=info)
140 6 : if (info /= 0) call judft_error("Can't allocate coul%mt2_c")
141 : endif
142 :
143 12 : if (.not. allocated(coul%mt3_c)) then
144 102 : allocate (coul%mt3_c(maxval(num_radbasfn) - 1, fi%atoms%nat, fi%atoms%nat), stat=info)
145 6 : if (info /= 0) call judft_error("Can't allocate coul%mt3_c")
146 : endif
147 12 : if (.not. coul%mtir%allocated()) call coul%mtir%alloc(.False., isize, isize)
148 : endif
149 36 : end subroutine t_coul_alloc
150 :
151 0 : subroutine t_coul_mini_alloc(coul, fi)
152 : implicit NONE
153 : type(t_fleurinput), intent(in) :: fi
154 : class(t_coul), intent(inout) :: coul
155 :
156 0 : if (.not. allocated(coul%mt1_r)) then
157 0 : allocate (coul%mt1_r(1,1,1,1))
158 : endif
159 0 : if (.not. allocated(coul%mt2_r)) then
160 0 : allocate (coul%mt2_r(1, 1, 1, 1))
161 : endif
162 :
163 0 : if (.not. allocated(coul%mt3_r)) then
164 0 : allocate (coul%mt3_r(1, 1, 1))
165 : endif
166 0 : if (.not. allocated(coul%mt1_c)) then
167 0 : allocate (coul%mt1_c(1, 1, 1, 1))
168 : endif
169 0 : if (.not. allocated(coul%mt2_c)) then
170 0 : allocate (coul%mt2_c(1, 1, 1, 1))
171 : endif
172 :
173 0 : if (.not. allocated(coul%mt3_c)) then
174 0 : allocate (coul%mt3_c(1, 1, 1))
175 : endif
176 :
177 0 : call coul%mtir%alloc(fi%sym%invs, 1, 1)
178 0 : end subroutine t_coul_mini_alloc
179 0 : end module m_types_coul
|