Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Auxiliary tools to redistribute cp_fm_type matrices before and after diagonalization.
10 : !> Heuristics are used to determine the optimal number of CPUs for diagonalization and the
11 : !> input matrices are redistributed if necessary
12 : !> \par History
13 : !> - [01.2018] moved redistribution related code from cp_fm_syevd here
14 : !> \author Nico Holmberg [01.2018]
15 : ! **************************************************************************************************
16 : MODULE cp_fm_diag_utils
17 : USE cp_blacs_env, ONLY: cp_blacs_env_create,&
18 : cp_blacs_env_release,&
19 : cp_blacs_env_type
20 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
21 : cp_fm_struct_release,&
22 : cp_fm_struct_type
23 : USE cp_fm_types, ONLY: cp_fm_create,&
24 : cp_fm_get_info,&
25 : cp_fm_release,&
26 : cp_fm_type
27 : USE cp_log_handling, ONLY: cp_get_default_logger,&
28 : cp_logger_get_default_io_unit,&
29 : cp_logger_type
30 : USE kinds, ONLY: dp
31 : USE mathlib, ONLY: gcd
32 : USE message_passing, ONLY: mp_para_env_type
33 : #include "../base/base_uses.f90"
34 :
35 : IMPLICIT NONE
36 :
37 : PRIVATE
38 :
39 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_fm_diag_utils'
40 :
41 : ! Information on redistribution
42 : TYPE, PUBLIC :: cp_fm_redistribute_info
43 : INTEGER :: matrix_order = -1
44 : INTEGER :: num_pe_old = -1 ! number of processes before a potential redistribute
45 : INTEGER :: num_pe_new = -1 ! number of processes after a potential redistribute
46 : INTEGER :: num_pe_opt = -1 ! optimal number of processes for the given matrix
47 : INTEGER :: num_pe_max_nz_col = -1 ! the maximal number of processes s.t. no column has zero width, may be < 0 if ignored
48 : LOGICAL :: redistribute = .FALSE. ! whether or not the matrix was actually redistributed
49 : CONTAINS
50 : PROCEDURE, PASS(self) :: write => cp_fm_redistribute_info_write
51 : END TYPE cp_fm_redistribute_info
52 :
53 : ! Container for redistribution settings and temporary work structs
54 : TYPE cp_fm_redistribute_type
55 : ! Settings
56 : INTEGER :: a = -1, x = -1
57 : LOGICAL :: should_print = .FALSE.
58 : LOGICAL :: elpa_force_redistribute = .FALSE.
59 : ! Temporaries
60 : INTEGER, DIMENSION(:), POINTER :: group_distribution => NULL(), &
61 : group_partition => NULL()
62 : TYPE(cp_blacs_env_type), POINTER :: blacs_env_new => NULL()
63 : TYPE(mp_para_env_type), POINTER :: para_env_new => NULL()
64 : END TYPE cp_fm_redistribute_type
65 :
66 : ! Permanent instance of the redistribute type
67 : TYPE(cp_fm_redistribute_type), PRIVATE, &
68 : SAVE :: work_redistribute
69 :
70 : ! Public subroutines
71 :
72 : PUBLIC :: cp_fm_redistribute_start, &
73 : cp_fm_redistribute_end, &
74 : cp_fm_redistribute_init
75 :
76 : CONTAINS
77 :
78 : ! **************************************************************************************************
79 : !> \brief Write the redistribute info nicely formatted to the given I/O unit
80 : !> \param self reference to the cp_fm_redistribute_info instance
81 : !> \param io_unit I/O unit to use for writing
82 : ! **************************************************************************************************
83 0 : SUBROUTINE cp_fm_redistribute_info_write(self, io_unit)
84 : CLASS(cp_fm_redistribute_info), INTENT(IN) :: self
85 : INTEGER, INTENT(IN) :: io_unit
86 :
87 0 : WRITE (UNIT=io_unit, FMT="(A)") ""
88 : WRITE (UNIT=io_unit, FMT="(T2,A,T71,I10)") &
89 0 : "CP_FM_DIAG| Number of processes over which the matrix is distributed ", self%num_pe_old, &
90 0 : "CP_FM_DIAG| Matrix order ", self%matrix_order
91 : WRITE (UNIT=io_unit, FMT="(T2,A,T71,I10)") &
92 0 : "CP_FM_DIAG| Optimal number of CPUs ", self%num_pe_opt
93 0 : IF (self%num_pe_max_nz_col < 0) THEN
94 : WRITE (UNIT=io_unit, FMT="(T2,A,T71,A10)") &
95 0 : "CP_FM_DIAG| Maximum number of CPUs (with non-zero columns) ", "<N/A>"
96 : ELSE
97 : WRITE (UNIT=io_unit, FMT="(T2,A,T71,I10)") &
98 0 : "CP_FM_DIAG| Maximum number of CPUs (with non-zero columns): ", self%num_pe_max_nz_col
99 : END IF
100 0 : IF (self%redistribute) THEN
101 : WRITE (UNIT=io_unit, FMT="(T2,A,T71,I10)") &
102 0 : "CP_FM_DIAG| Number of processes for the redistribution ", self%num_pe_new
103 : ELSE
104 : WRITE (UNIT=io_unit, FMT="(T2,A)") &
105 0 : "CP_FM_DIAG| The matrix will NOT be redistributed"
106 : END IF
107 0 : WRITE (UNIT=io_unit, FMT="(A)") ""
108 :
109 0 : END SUBROUTINE cp_fm_redistribute_info_write
110 :
111 : ! **************************************************************************************************
112 : !> \brief Releases the temporary storage needed when redistributing arrays
113 : !> \param has_redistributed flag that determines if the processors holds a part of the
114 : !> redistributed array
115 : !> \author Nico Holmberg [01.2018]
116 : ! **************************************************************************************************
117 268841 : SUBROUTINE cp_fm_redistribute_work_finalize(has_redistributed)
118 : LOGICAL, INTENT(IN) :: has_redistributed
119 :
120 268841 : IF (ASSOCIATED(work_redistribute%group_distribution)) THEN
121 268841 : IF (has_redistributed) THEN
122 136294 : CALL cp_blacs_env_release(work_redistribute%blacs_env_new)
123 : END IF
124 268841 : CALL work_redistribute%para_env_new%free()
125 268841 : DEALLOCATE (work_redistribute%para_env_new)
126 268841 : DEALLOCATE (work_redistribute%group_distribution)
127 268841 : DEALLOCATE (work_redistribute%group_partition)
128 : END IF
129 : ! Return work to its initial state
130 268841 : work_redistribute = cp_fm_redistribute_type()
131 :
132 268841 : END SUBROUTINE cp_fm_redistribute_work_finalize
133 :
134 : ! **************************************************************************************************
135 : !> \brief Initializes the parameters that determine how to calculate the optimal number of CPUs
136 : !> for diagonalizing a matrix. The parameters are read from the GLOBAL input section.
137 : !> \param a integer parameter used to define the rule for determining the optimal
138 : !> number of CPUs for diagonalization
139 : !> \param x integer parameter used to define the rule for determining the optimal
140 : !> number of CPUs for diagonalization
141 : !> \param should_print flag that determines if information about the redistribution process
142 : !> should be printed
143 : !> \param elpa_force_redistribute flag that if redistribution should always be performed when
144 : !> the ELPA diagonalization library is in use
145 : !> \author Nico Holmberg [01.2018]
146 : ! **************************************************************************************************
147 11087 : SUBROUTINE cp_fm_redistribute_init(a, x, should_print, elpa_force_redistribute)
148 : INTEGER, INTENT(IN) :: a, x
149 : LOGICAL, INTENT(IN) :: should_print, elpa_force_redistribute
150 :
151 : work_redistribute%a = a
152 : work_redistribute%x = x
153 : work_redistribute%should_print = should_print
154 : work_redistribute%elpa_force_redistribute = elpa_force_redistribute
155 : ! Init work
156 11087 : work_redistribute = cp_fm_redistribute_type()
157 :
158 11087 : END SUBROUTINE cp_fm_redistribute_init
159 :
160 : ! **************************************************************************************************
161 : !> \brief Calculates the optimal number of CPUs for diagonalizing a matrix.
162 : !> \param size the size of the diagonalized matrix
163 : !> \return the optimal number of CPUs
164 : !> \author Nico Holmberg [01.2018]
165 : ! **************************************************************************************************
166 272985 : PURE FUNCTION cp_fm_diag_get_optimal_ncpu(size) RESULT(ncpu)
167 : INTEGER, INTENT(IN) :: size
168 : INTEGER :: ncpu
169 :
170 : ncpu = ((size + work_redistribute%a*work_redistribute%x - 1)/ &
171 272985 : (work_redistribute%a*work_redistribute%x))*work_redistribute%a
172 :
173 272985 : END FUNCTION cp_fm_diag_get_optimal_ncpu
174 :
175 : #if defined(__parallel)
176 : ! **************************************************************************************************
177 : !> \brief Determines the largest number of CPUs a matrix can be distributed on without any of the
178 : !> processors getting a zero-width column (currently only needed for ELPA).
179 : !> \param matrix the matrix that will be diagonalized
180 : !> \return the maximum number of CPUs for ELPA
181 : !> \author Nico Holmberg [01.2018]
182 : ! **************************************************************************************************
183 16222 : FUNCTION cp_fm_max_ncpu_non_zero_column(matrix) RESULT(ncpu)
184 : TYPE(cp_fm_type), INTENT(IN) :: matrix
185 : INTEGER :: ncpu
186 :
187 : INTEGER :: gcd_max, ipe, jpe, ncol_block, &
188 : ncol_global, npcol, nrow_block, &
189 : nrow_global, num_pe_old, nzero
190 16222 : INTEGER, DIMENSION(:), POINTER :: ncol_locals
191 : INTEGER, EXTERNAL :: numroc
192 :
193 16222 : NULLIFY (ncol_locals)
194 : ! First check if there are any zero width columns in current layout
195 : CALL cp_fm_get_info(matrix, ncol_locals=ncol_locals, &
196 : nrow_global=nrow_global, ncol_global=ncol_global, &
197 16222 : nrow_block=nrow_block, ncol_block=ncol_block)
198 32444 : nzero = COUNT(ncol_locals == 0)
199 16222 : num_pe_old = matrix%matrix_struct%para_env%num_pe
200 16222 : ncpu = num_pe_old - nzero
201 :
202 : ! Avoid layouts with odd number of CPUs (blacs grid layout will be square)
203 16222 : IF (ncpu > 2) THEN
204 0 : ncpu = ncpu - MODULO(ncpu, 2)
205 : END IF
206 :
207 : ! if there are no zero-width columns and the number of processors was even, leave it at that
208 16222 : IF (ncpu == num_pe_old) THEN
209 : RETURN
210 : END IF
211 :
212 : ! Iteratively search for the maximum number of CPUs for ELPA
213 : ! On each step, we test whether the blacs grid created with ncpu processes
214 : ! contains any columns with zero width
215 0 : DO WHILE (ncpu > 1)
216 : ! Determine layout of new blacs grid with ncpu CPUs
217 : ! (snippet copied from cp_blacs_env.F:cp_blacs_env_create)
218 0 : gcd_max = -1
219 0 : DO ipe = 1, CEILING(SQRT(REAL(ncpu, dp)))
220 0 : jpe = ncpu/ipe
221 0 : IF (ipe*jpe /= ncpu) THEN
222 : CYCLE
223 : END IF
224 0 : IF (gcd(ipe, jpe) >= gcd_max) THEN
225 0 : npcol = jpe
226 0 : gcd_max = gcd(ipe, jpe)
227 : END IF
228 : END DO
229 :
230 : ! Count the number of processors without any columns
231 : ! (snippet copied from cp_fm_struct.F:cp_fm_struct_create)
232 0 : nzero = 0
233 0 : DO ipe = 0, npcol - 1
234 0 : IF (numroc(ncol_global, ncol_block, ipe, 0, npcol) == 0) THEN
235 0 : nzero = nzero + 1
236 : END IF
237 : END DO
238 :
239 0 : IF (nzero == 0) THEN
240 : EXIT
241 : END IF
242 :
243 0 : ncpu = ncpu - nzero
244 :
245 0 : IF (ncpu > 2) THEN
246 0 : ncpu = ncpu - MODULO(ncpu, 2)
247 : END IF
248 : END DO
249 :
250 16222 : END FUNCTION cp_fm_max_ncpu_non_zero_column
251 : #endif
252 :
253 : ! **************************************************************************************************
254 : !> \brief Determines the optimal number of CPUs for matrix diagonalization and redistributes
255 : !> the input matrices if necessary
256 : !> \param matrix the input cp_fm_type matrix to be diagonalized
257 : !> \param eigenvectors the cp_fm_type matrix that will hold the eigenvectors of the input matrix
258 : !> \param matrix_new the redistributed input matrix which will subsequently be diagonalized,
259 : !> or a pointer to the original matrix if no redistribution is required
260 : !> \param eigenvectors_new the redistributed eigenvectors matrix, or a pointer to the original
261 : !> matrix if no redistribution is required
262 : !> \param caller_is_elpa flag that determines if ELPA is used for diagonalization
263 : !> \param redist_info get info about the redistribution
264 : !> \par History
265 : !> - [01.2018] created by moving redistribution related code from cp_fm_syevd here
266 : !> \author Nico Holmberg [01.2018]
267 : ! **************************************************************************************************
268 16222 : SUBROUTINE cp_fm_redistribute_start(matrix, eigenvectors, matrix_new, eigenvectors_new, &
269 : caller_is_elpa, redist_info)
270 :
271 : TYPE(cp_fm_type), INTENT(IN) :: matrix, eigenvectors
272 : TYPE(cp_fm_type), INTENT(OUT) :: matrix_new, eigenvectors_new
273 : LOGICAL, OPTIONAL, INTENT(IN) :: caller_is_elpa
274 :
275 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_redistribute_start'
276 :
277 : INTEGER :: handle
278 : LOGICAL :: is_elpa
279 : TYPE(cp_fm_redistribute_info), OPTIONAL, INTENT(OUT) :: redist_info
280 :
281 : #if defined(__parallel)
282 : REAL(KIND=dp) :: fake_local_data(1, 1)
283 : INTEGER :: fake_descriptor(9), mepos_old, &
284 : io_unit, ngroups, ncol_block, blksize, nrow_block
285 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_new
286 : TYPE(mp_para_env_type), POINTER :: para_env
287 : TYPE(cp_logger_type), POINTER :: logger
288 : TYPE(cp_fm_redistribute_info) :: rdinfo
289 : #endif
290 :
291 272985 : CALL timeset(routineN, handle)
292 272985 : is_elpa = .FALSE.
293 272985 : IF (PRESENT(caller_is_elpa)) THEN
294 : #if defined(__ELPA)
295 16222 : is_elpa = caller_is_elpa
296 : #else
297 : CPABORT("CP2K compiled without the ELPA library.")
298 : #endif
299 : END IF
300 :
301 : #if defined(__parallel)
302 :
303 272985 : logger => cp_get_default_logger()
304 272985 : io_unit = cp_logger_get_default_io_unit(logger)
305 :
306 : ! first figure out the optimal number of cpus
307 : ! this is pure heuristics, the defaults are based on rosa timings
308 : ! that demonstrate that timings go up sharply if too many tasks are used
309 : ! we take a multiple of 4, and approximately n/60
310 272985 : para_env => matrix%matrix_struct%para_env
311 272985 : mepos_old = para_env%mepos
312 272985 : ncol_block = -1 ! normally we also want to adjust the block size according to the optimal # of CPUs
313 272985 : nrow_block = -1
314 272985 : blksize = -1
315 :
316 272985 : rdinfo%matrix_order = matrix%matrix_struct%nrow_global
317 272985 : rdinfo%num_pe_old = para_env%num_pe
318 272985 : rdinfo%num_pe_opt = cp_fm_diag_get_optimal_ncpu(rdinfo%matrix_order)
319 272985 : rdinfo%num_pe_new = rdinfo%num_pe_opt
320 : rdinfo%num_pe_max_nz_col = -1
321 : rdinfo%redistribute = .FALSE.
322 :
323 272985 : IF (is_elpa) THEN
324 : ! with ELPA we don't have to redistribute if not necessary (scales, unlike ScaLAPACK)
325 16222 : rdinfo%num_pe_new = rdinfo%num_pe_old
326 :
327 : ! BUT: Diagonalization with ELPA fails when a processor column has zero width
328 : ! Determine the maximum number of CPUs the matrix can be distributed without zero-width columns
329 : ! for the current block size.
330 16222 : rdinfo%num_pe_max_nz_col = cp_fm_max_ncpu_non_zero_column(matrix)
331 :
332 : ! if the user wants to redistribute to the ScaLAPACK optimal number of CPUs anyway, let him if it's safe.
333 16222 : IF (work_redistribute%elpa_force_redistribute .AND. rdinfo%num_pe_opt < rdinfo%num_pe_max_nz_col) THEN
334 : ! Use heuristics to determine the need for redistribution (when num_pe_opt is smaller than the safe maximum)
335 : ! in this case we can also take the block size used for ScaLAPACK
336 0 : rdinfo%num_pe_new = rdinfo%num_pe_opt
337 16222 : ELSE IF (rdinfo%num_pe_old > rdinfo%num_pe_max_nz_col) THEN
338 : ! Otherwise, only redistribute if we have to
339 0 : rdinfo%num_pe_new = rdinfo%num_pe_max_nz_col
340 : ! do NOT let cp_fm_struct_create automatically adjust the block size because the
341 : ! calculated number of processors such that no block has 0 columns wouldn't match (see #578):
342 : ! if the automatically chosen block size is larger than the present one we would still end
343 : ! up with empty processors
344 : END IF
345 :
346 16222 : CALL cp_fm_get_info(matrix, ncol_block=ncol_block, nrow_block=nrow_block)
347 :
348 : ! On GPUs, ELPA requires the block size to be a power of 2
349 16222 : blksize = 1
350 99498 : DO WHILE (2*blksize <= MIN(nrow_block, ncol_block))
351 16222 : blksize = blksize*2
352 : END DO
353 16222 : nrow_block = blksize
354 16222 : ncol_block = blksize
355 : END IF
356 :
357 : ! finally, only redistribute if we're going to use less CPUs than before or changed the block size
358 : rdinfo%redistribute = (rdinfo%num_pe_old > rdinfo%num_pe_new) .OR. (blksize >= 0 .AND. &
359 272985 : ((blksize /= matrix%matrix_struct%ncol_block) .OR. (blksize /= matrix%matrix_struct%nrow_block)))
360 :
361 272985 : IF (work_redistribute%should_print .AND. io_unit > 0) THEN
362 0 : IF (is_elpa) THEN
363 0 : IF (work_redistribute%elpa_force_redistribute) THEN
364 : WRITE (UNIT=io_unit, FMT="(T2,A,T78,A3)") &
365 0 : "CP_FM_DIAG| Force redistribute (ELPA):", "YES"
366 : ELSE
367 : WRITE (UNIT=io_unit, FMT="(T2,A,T79,A2)") &
368 0 : "CP_FM_DIAG| Force redistribute (ELPA):", "NO"
369 : END IF
370 : END IF
371 0 : CALL rdinfo%write(io_unit)
372 : END IF
373 272985 : CALL para_env%sync()
374 :
375 : ! if the optimal is smaller than num_pe, we will redistribute the input matrix
376 272985 : IF (rdinfo%redistribute) THEN
377 : ! split comm, the first num_pe_new tasks will do the work
378 806523 : ALLOCATE (work_redistribute%group_distribution(0:rdinfo%num_pe_old - 1))
379 268841 : ALLOCATE (work_redistribute%group_partition(0:1))
380 806523 : work_redistribute%group_partition = [rdinfo%num_pe_new, rdinfo%num_pe_old - rdinfo%num_pe_new]
381 268841 : ALLOCATE (work_redistribute%para_env_new)
382 : CALL work_redistribute%para_env_new%from_split( &
383 : comm=para_env, ngroups=ngroups, group_distribution=work_redistribute%group_distribution, &
384 268841 : n_subgroups=2, group_partition=work_redistribute%group_partition)
385 :
386 268841 : IF (work_redistribute%group_distribution(mepos_old) == 0) THEN
387 :
388 : ! create blacs, should inherit the preferences for the layout and so on, from the higher level
389 136294 : NULLIFY (work_redistribute%blacs_env_new)
390 136294 : CALL cp_blacs_env_create(blacs_env=work_redistribute%blacs_env_new, para_env=work_redistribute%para_env_new)
391 :
392 : ! create new matrix
393 136294 : NULLIFY (fm_struct_new)
394 136294 : IF (nrow_block == -1 .OR. ncol_block == -1) THEN
395 : CALL cp_fm_struct_create(fmstruct=fm_struct_new, &
396 : para_env=work_redistribute%para_env_new, &
397 : context=work_redistribute%blacs_env_new, &
398 : nrow_global=rdinfo%matrix_order, ncol_global=rdinfo%matrix_order, &
399 130255 : ncol_block=ncol_block, nrow_block=nrow_block)
400 : ELSE
401 : CALL cp_fm_struct_create(fmstruct=fm_struct_new, &
402 : para_env=work_redistribute%para_env_new, &
403 : context=work_redistribute%blacs_env_new, &
404 : nrow_global=rdinfo%matrix_order, ncol_global=rdinfo%matrix_order, &
405 6039 : ncol_block=ncol_block, nrow_block=nrow_block, force_block=.TRUE.)
406 : END IF
407 136294 : CALL cp_fm_create(matrix_new, matrix_struct=fm_struct_new, name="yevd_new_mat")
408 136294 : CALL cp_fm_create(eigenvectors_new, matrix_struct=fm_struct_new, name="yevd_new_vec")
409 136294 : CALL cp_fm_struct_release(fm_struct_new)
410 :
411 : ! redistribute old
412 : CALL pdgemr2d(rdinfo%matrix_order, rdinfo%matrix_order, matrix%local_data(1, 1), 1, 1, &
413 : matrix%matrix_struct%descriptor, &
414 : matrix_new%local_data(1, 1), 1, 1, matrix_new%matrix_struct%descriptor, &
415 136294 : matrix%matrix_struct%context)
416 : ELSE
417 : ! these tasks must help redistribute (they own part of the data),
418 : ! but need fake 'new' data, and their descriptor must indicate this with -1
419 : ! see also scalapack comments on pdgemr2d
420 1325470 : fake_descriptor = -1
421 : CALL pdgemr2d(rdinfo%matrix_order, rdinfo%matrix_order, matrix%local_data(1, 1), 1, 1, &
422 : matrix%matrix_struct%descriptor, &
423 : fake_local_data(1, 1), 1, 1, fake_descriptor, &
424 132547 : matrix%matrix_struct%context)
425 : END IF
426 : ELSE
427 : ! No need to redistribute, just return pointers to the original arrays
428 4144 : matrix_new = matrix
429 4144 : eigenvectors_new = eigenvectors
430 : END IF
431 :
432 272985 : IF (PRESENT(redist_info)) THEN
433 16222 : redist_info = rdinfo
434 : END IF
435 : #else
436 :
437 : MARK_USED(matrix)
438 : MARK_USED(eigenvectors)
439 : MARK_USED(matrix_new)
440 : MARK_USED(eigenvectors_new)
441 : MARK_USED(redist_info)
442 : CPABORT("Routine called in non-parallel case.")
443 : #endif
444 :
445 272985 : CALL timestop(handle)
446 :
447 272985 : END SUBROUTINE cp_fm_redistribute_start
448 :
449 : ! **************************************************************************************************
450 : !> \brief Redistributes eigenvectors and eigenvalues back to the original communicator group
451 : !> \param matrix the input cp_fm_type matrix to be diagonalized
452 : !> \param eigenvectors the cp_fm_type matrix that will hold the eigenvectors of the input matrix
453 : !> \param eig global array holding the eigenvalues of the input matrixmatrix
454 : !> \param matrix_new the redistributed input matrix which will subsequently be diagonalized,
455 : !> or a pointer to the original matrix if no redistribution is required
456 : !> \param eigenvectors_new the redistributed eigenvectors matrix, or a pointer to the original
457 : !> matrix if no redistribution is required
458 : !> \par History
459 : !> - [01.2018] created by moving redistribution related code from cp_fm_syevd here
460 : !> \author Nico Holmberg [01.2018]
461 : ! **************************************************************************************************
462 272985 : SUBROUTINE cp_fm_redistribute_end(matrix, eigenvectors, eig, matrix_new, eigenvectors_new)
463 :
464 : TYPE(cp_fm_type), INTENT(IN) :: matrix, eigenvectors
465 : REAL(KIND=dp), DIMENSION(:), INTENT(INOUT) :: eig
466 : TYPE(cp_fm_type), INTENT(INOUT) :: matrix_new, eigenvectors_new
467 :
468 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_redistribute_end'
469 :
470 : INTEGER :: handle
471 : #if defined(__parallel)
472 : REAL(KIND=dp) :: fake_local_data(1, 1)
473 : INTEGER :: fake_descriptor(9), mepos_old, n
474 : TYPE(mp_para_env_type), POINTER :: para_env
475 : #endif
476 :
477 272985 : CALL timeset(routineN, handle)
478 :
479 : #if defined(__parallel)
480 :
481 : ! Check if matrix was redistributed
482 272985 : IF (ASSOCIATED(work_redistribute%group_distribution)) THEN
483 268841 : n = matrix%matrix_struct%nrow_global
484 268841 : para_env => matrix%matrix_struct%para_env
485 268841 : mepos_old = para_env%mepos
486 :
487 268841 : IF (work_redistribute%group_distribution(mepos_old) == 0) THEN
488 : ! redistribute results on CPUs that hold the redistributed matrix
489 : CALL pdgemr2d(n, n, eigenvectors_new%local_data(1, 1), 1, 1, eigenvectors_new%matrix_struct%descriptor, &
490 : eigenvectors%local_data(1, 1), 1, 1, eigenvectors%matrix_struct%descriptor, &
491 136294 : eigenvectors%matrix_struct%context)
492 136294 : CALL cp_fm_release(matrix_new)
493 136294 : CALL cp_fm_release(eigenvectors_new)
494 : ELSE
495 : ! these tasks must help redistribute (they own part of the data),
496 : ! but need fake 'new' data, and their descriptor must indicate this with -1
497 : ! see also scalapack comments on pdgemr2d
498 1325470 : fake_descriptor = -1
499 : CALL pdgemr2d(n, n, fake_local_data(1, 1), 1, 1, fake_descriptor, &
500 : eigenvectors%local_data(1, 1), 1, 1, eigenvectors%matrix_struct%descriptor, &
501 132547 : eigenvectors%matrix_struct%context)
502 : END IF
503 : ! free work
504 268841 : CALL cp_fm_redistribute_work_finalize(work_redistribute%group_distribution(mepos_old) == 0)
505 :
506 : ! finally, also the eigenvalues need to end up on the non-group member tasks
507 7710393 : CALL para_env%bcast(eig, 0)
508 : END IF
509 :
510 : #else
511 :
512 : MARK_USED(matrix)
513 : MARK_USED(eigenvectors)
514 : MARK_USED(eig)
515 : MARK_USED(matrix_new)
516 : MARK_USED(eigenvectors_new)
517 : CPABORT("Routine called in non-parallel case.")
518 : #endif
519 :
520 272985 : CALL timestop(handle)
521 :
522 272985 : END SUBROUTINE cp_fm_redistribute_end
523 :
524 0 : END MODULE cp_fm_diag_utils
|
