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 DBCSR operations in CP2K
10 : !> \author Urban Borstnik
11 : !> \date 2009-05-12
12 : !> \version 0.8
13 : !>
14 : !> <b>Modification history:</b>
15 : !> - Created 2009-05-12
16 : !> - Generalized sm_fm_mulitply for matrices w/ different row/col block size (A. Bussy, 11.2018)
17 : ! **************************************************************************************************
18 : MODULE cp_dbcsr_operations
19 : USE cp_blacs_env, ONLY: cp_blacs_env_type
20 : USE cp_dbcsr_api, ONLY: &
21 : dbcsr_add, dbcsr_complete_redistribute, dbcsr_convert_sizes_to_offsets, dbcsr_copy, &
22 : dbcsr_create, dbcsr_deallocate_matrix, dbcsr_desymmetrize, dbcsr_distribution_get, &
23 : dbcsr_distribution_new, dbcsr_distribution_release, dbcsr_distribution_type, &
24 : dbcsr_get_info, dbcsr_get_matrix_type, dbcsr_iterator_blocks_left, &
25 : dbcsr_iterator_next_block, dbcsr_iterator_readonly_start, dbcsr_iterator_start, &
26 : dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, dbcsr_p_type, dbcsr_release, &
27 : dbcsr_scale, dbcsr_type, dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, &
28 : dbcsr_type_symmetric, dbcsr_valid_index, dbcsr_verify_matrix
29 : USE cp_dbcsr_contrib, ONLY: dbcsr_frobenius_norm,&
30 : dbcsr_reserve_all_blocks
31 : USE cp_fm_basic_linalg, ONLY: cp_fm_gemm
32 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
33 : cp_fm_struct_release,&
34 : cp_fm_struct_type
35 : USE cp_fm_types, ONLY: cp_fm_create,&
36 : cp_fm_get_info,&
37 : cp_fm_release,&
38 : cp_fm_to_fm,&
39 : cp_fm_type
40 : USE distribution_2d_types, ONLY: distribution_2d_get,&
41 : distribution_2d_type
42 : USE kinds, ONLY: default_string_length,&
43 : dp,&
44 : int_8
45 : USE mathlib, ONLY: gcd,&
46 : lcm
47 : USE message_passing, ONLY: mp_para_env_type
48 :
49 : !$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num, omp_get_num_threads
50 : #include "base/base_uses.f90"
51 :
52 : IMPLICIT NONE
53 : PRIVATE
54 :
55 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_operations'
56 : LOGICAL, PARAMETER :: debug_mod = .FALSE.
57 :
58 : INTEGER, SAVE, PUBLIC :: max_elements_per_block = 32
59 :
60 : ! A plan belongs to a fixed S/A workspace and FM distribution. Only values may
61 : ! change during its lifetime. Its owner must release it before the workspace or
62 : ! FM distribution changes. All owned storage is allocatable. No matrix data
63 : ! pointers or MPI requests are retained. Calls using one plan must be serialized.
64 : TYPE dbcsr_fm_thread_plan_type
65 : INTEGER, ALLOCATABLE :: blocks(:, :), slots(:)
66 : END TYPE dbcsr_fm_thread_plan_type
67 :
68 : TYPE, PUBLIC :: dbcsr_fm_copy_plan_type
69 : PRIVATE
70 : TYPE(cp_fm_struct_type), POINTER :: matrix_struct => NULL()
71 : CHARACTER :: matrix_type = " "
72 : TYPE(dbcsr_fm_thread_plan_type), ALLOCATABLE :: thread(:)
73 : INTEGER, ALLOCATABLE :: send_count(:), send_disp(:), recv_count(:), recv_disp(:), recv_index(:, :)
74 : REAL(KIND=dp), ALLOCATABLE :: send_values(:), recv_values(:)
75 : END TYPE dbcsr_fm_copy_plan_type
76 :
77 : PUBLIC :: dbcsr_multiply_local
78 :
79 : ! CP2K API emulation
80 : PUBLIC :: copy_fm_to_dbcsr, copy_dbcsr_to_fm, &
81 : cp_dbcsr_sm_fm_multiply, cp_dbcsr_plus_fm_fm_t, &
82 : copy_dbcsr_to_fm_bc, copy_fm_to_dbcsr_bc, cp_fm_to_dbcsr_row_template, &
83 : cp_dbcsr_m_by_n_from_template, cp_dbcsr_m_by_n_from_row_template, &
84 : dbcsr_create_dist_r_unrot
85 :
86 : ! distribution_2d_type compatibility
87 : PUBLIC :: cp_dbcsr_dist2d_to_dist
88 :
89 : PUBLIC :: dbcsr_copy_columns_hack
90 :
91 : ! matrix set
92 : PUBLIC :: dbcsr_allocate_matrix_set
93 : PUBLIC :: dbcsr_deallocate_matrix_set
94 :
95 : INTERFACE dbcsr_allocate_matrix_set
96 : MODULE PROCEDURE allocate_dbcsr_matrix_set_1d
97 : MODULE PROCEDURE allocate_dbcsr_matrix_set_2d
98 : MODULE PROCEDURE allocate_dbcsr_matrix_set_3d
99 : MODULE PROCEDURE allocate_dbcsr_matrix_set_4d
100 : MODULE PROCEDURE allocate_dbcsr_matrix_set_5d
101 : END INTERFACE
102 :
103 : INTERFACE dbcsr_deallocate_matrix_set
104 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_1d
105 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_2d
106 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_3d
107 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_4d
108 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_5d
109 : END INTERFACE
110 :
111 : CONTAINS
112 :
113 : ! **************************************************************************************************
114 : !> \brief Copy a BLACS matrix to a dbcsr matrix.
115 : !>
116 : !> real_matrix=beta*real_matrix+alpha*fm
117 : !> beta defaults to 0, alpha to 1
118 : !> \param[in] fm full matrix
119 : !> \param[out] matrix DBCSR matrix
120 : !> \param[in] keep_sparsity (optional) retains the sparsity of the input
121 : !> matrix
122 : !> \date 2009-10-13
123 : !> \par History
124 : !> 2009-10-13 rewritten based on copy_dbcsr_to_fm
125 : !> \author Urban Borstnik
126 : !> \version 2.0
127 : ! **************************************************************************************************
128 1623818 : SUBROUTINE copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
129 : TYPE(cp_fm_type), INTENT(IN) :: fm
130 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix
131 : LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
132 :
133 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_fm_to_dbcsr'
134 :
135 : INTEGER :: handle
136 : LOGICAL :: my_keep_sparsity
137 : TYPE(dbcsr_type) :: bc_mat, redist_mat
138 :
139 1623818 : CALL timeset(routineN, handle)
140 :
141 1623818 : my_keep_sparsity = .FALSE.
142 1623818 : IF (PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
143 :
144 1623818 : CALL copy_fm_to_dbcsr_bc(fm, bc_mat)
145 :
146 1623818 : IF (my_keep_sparsity) THEN
147 296140 : CALL dbcsr_create(redist_mat, template=matrix)
148 296140 : CALL dbcsr_complete_redistribute(bc_mat, redist_mat)
149 296140 : CALL dbcsr_copy(matrix, redist_mat, keep_sparsity=.TRUE.)
150 296140 : CALL dbcsr_release(redist_mat)
151 : ELSE
152 1327678 : CALL dbcsr_complete_redistribute(bc_mat, matrix)
153 : END IF
154 :
155 1623818 : CALL dbcsr_release(bc_mat)
156 :
157 1623818 : CALL timestop(handle)
158 1623818 : END SUBROUTINE copy_fm_to_dbcsr
159 :
160 : ! **************************************************************************************************
161 : !> \brief Copy a BLACS matrix to a dbcsr matrix with a special block-cyclic distribution,
162 : !> which requires no complete redistribution.
163 : !> \param fm ...
164 : !> \param bc_mat ...
165 : ! **************************************************************************************************
166 1631060 : SUBROUTINE copy_fm_to_dbcsr_bc(fm, bc_mat)
167 : TYPE(cp_fm_type), INTENT(IN) :: fm
168 : TYPE(dbcsr_type), INTENT(INOUT) :: bc_mat
169 :
170 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_fm_to_dbcsr_bc'
171 :
172 : INTEGER :: col, handle, ncol_block, ncol_global, &
173 : nrow_block, nrow_global, row
174 1631060 : INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
175 1631060 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
176 1631060 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
177 1631060 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: dbcsr_block, fm_block
178 : TYPE(dbcsr_distribution_type) :: bc_dist
179 : TYPE(dbcsr_iterator_type) :: iter
180 :
181 1631060 : CALL timeset(routineN, handle)
182 :
183 : ! Create processor grid
184 1631060 : pgrid => fm%matrix_struct%context%blacs2mpi
185 :
186 : ! Create a block-cyclic distribution compatible with the FM matrix.
187 1631060 : nrow_block = fm%matrix_struct%nrow_block
188 1631060 : ncol_block = fm%matrix_struct%ncol_block
189 1631060 : nrow_global = fm%matrix_struct%nrow_global
190 1631060 : ncol_global = fm%matrix_struct%ncol_global
191 1631060 : NULLIFY (col_blk_size, row_blk_size)
192 : CALL dbcsr_create_dist_block_cyclic(bc_dist, &
193 : nrows=nrow_global, ncolumns=ncol_global, & ! Actual full matrix size
194 : nrow_block=nrow_block, ncol_block=ncol_block, & ! BLACS parameters
195 : group_handle=fm%matrix_struct%para_env%get_handle(), pgrid=pgrid, &
196 1631060 : row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size) ! block-cyclic row/col sizes
197 :
198 : ! Create the block-cyclic DBCSR matrix
199 : CALL dbcsr_create(bc_mat, "Block-cyclic ", bc_dist, &
200 1631060 : dbcsr_type_no_symmetry, row_blk_size, col_blk_size, reuse_arrays=.TRUE.)
201 1631060 : CALL dbcsr_distribution_release(bc_dist)
202 :
203 : ! allocate all blocks
204 1631060 : CALL dbcsr_reserve_all_blocks(bc_mat)
205 :
206 1631060 : CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
207 :
208 : ! Copy the FM data to the block-cyclic DBCSR matrix. This step
209 : ! could be skipped with appropriate DBCSR index manipulation.
210 1631060 : fm_block => fm%local_data
211 : !$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
212 1631060 : !$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
213 : CALL dbcsr_iterator_start(iter, bc_mat)
214 : DO WHILE (dbcsr_iterator_blocks_left(iter))
215 : CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
216 : dbcsr_block(:, :) = fm_block(first_row(row):last_row(row), first_col(col):last_col(col))
217 : END DO
218 : CALL dbcsr_iterator_stop(iter)
219 : !$OMP END PARALLEL
220 :
221 1631060 : CALL timestop(handle)
222 3262120 : END SUBROUTINE copy_fm_to_dbcsr_bc
223 :
224 : ! **************************************************************************************************
225 : !> \brief Copy a DBCSR matrix to a BLACS matrix
226 : !> \param[in] matrix DBCSR matrix
227 : !> \param[out] fm full matrix
228 : !> \param plan optional S/A routing plan; source block layout, FM structure and
229 : !> OpenMP thread distribution must stay fixed. Discard it before they change.
230 : ! **************************************************************************************************
231 1252401 : SUBROUTINE copy_dbcsr_to_fm(matrix, fm, plan)
232 : TYPE(dbcsr_type), INTENT(IN) :: matrix
233 : TYPE(cp_fm_type), INTENT(INOUT) :: fm
234 : TYPE(dbcsr_fm_copy_plan_type), INTENT(INOUT), &
235 : OPTIONAL :: plan
236 :
237 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_dbcsr_to_fm'
238 :
239 : CHARACTER(len=default_string_length) :: name
240 : INTEGER :: group_handle, handle, ncol_block, &
241 : nfullcols_total, nfullrows_total, &
242 : nrow_block
243 1252401 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
244 1252401 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
245 : TYPE(dbcsr_distribution_type) :: bc_dist, dist
246 1252401 : TYPE(dbcsr_fm_copy_plan_type) :: local_plan
247 : TYPE(dbcsr_type) :: bc_mat, matrix_nosym
248 :
249 1252401 : CALL timeset(routineN, handle)
250 :
251 : ! check compatibility
252 : CALL dbcsr_get_info(matrix, &
253 : name=name, &
254 : distribution=dist, &
255 : nfullrows_total=nfullrows_total, &
256 1252401 : nfullcols_total=nfullcols_total)
257 :
258 1252401 : CPASSERT(fm%matrix_struct%nrow_global == nfullrows_total)
259 1252401 : CPASSERT(fm%matrix_struct%ncol_global == nfullcols_total)
260 :
261 365401 : SELECT CASE (dbcsr_get_matrix_type(matrix))
262 : CASE (dbcsr_type_symmetric, dbcsr_type_antisymmetric)
263 365401 : IF (PRESENT(plan)) THEN
264 64 : CALL copy_dbcsr_to_fm_sym(matrix, fm, plan)
265 : ELSE
266 365337 : CALL copy_dbcsr_to_fm_sym(matrix, fm, local_plan)
267 : END IF
268 365401 : CALL timestop(handle)
269 1617802 : RETURN
270 : END SELECT
271 :
272 : ! Reusable plans are only defined for the S/A direct-copy path.
273 887000 : CPASSERT(.NOT. PRESENT(plan))
274 :
275 : ! info about the full matrix
276 887000 : nrow_block = fm%matrix_struct%nrow_block
277 887000 : ncol_block = fm%matrix_struct%ncol_block
278 :
279 : ! Convert DBCSR to a block-cyclic
280 887000 : NULLIFY (col_blk_size, row_blk_size)
281 887000 : CALL dbcsr_distribution_get(dist, group=group_handle, pgrid=pgrid)
282 : CALL dbcsr_create_dist_block_cyclic(bc_dist, &
283 : nrows=nfullrows_total, ncolumns=nfullcols_total, &
284 : nrow_block=nrow_block, ncol_block=ncol_block, &
285 : group_handle=group_handle, pgrid=pgrid, &
286 887000 : row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
287 :
288 : CALL dbcsr_create(bc_mat, "Block-cyclic"//name, bc_dist, &
289 887000 : dbcsr_type_no_symmetry, row_blk_size, col_blk_size, reuse_arrays=.TRUE.)
290 887000 : CALL dbcsr_distribution_release(bc_dist)
291 :
292 887000 : IF (dbcsr_get_matrix_type(matrix) == dbcsr_type_no_symmetry) THEN
293 887000 : CALL dbcsr_complete_redistribute(matrix, bc_mat)
294 : ELSE
295 0 : CALL dbcsr_create(matrix_nosym, template=matrix, matrix_type=dbcsr_type_no_symmetry)
296 0 : CALL dbcsr_desymmetrize(matrix, matrix_nosym)
297 0 : CALL dbcsr_complete_redistribute(matrix_nosym, bc_mat)
298 887000 : CALL dbcsr_release(matrix_nosym)
299 : END IF
300 :
301 887000 : CALL copy_dbcsr_to_fm_bc(bc_mat, fm)
302 :
303 887000 : CALL dbcsr_release(bc_mat)
304 :
305 887000 : CALL timestop(handle)
306 2870139 : END SUBROUTINE copy_dbcsr_to_fm
307 :
308 : ! **************************************************************************************************
309 : !> \brief Send S/A values using a fixed routing plan, built on first use.
310 : !> Validate the block sequence while packing, without adding a collective
311 : !> for cache validation. Structural changes require a fresh caller-owned plan.
312 : !> \param matrix symmetric or antisymmetric DBCSR workspace
313 : !> \param fm full matrix on the source communicator
314 : !> \param plan reusable routing and buffers; no outstanding communication on return
315 : ! **************************************************************************************************
316 365401 : SUBROUTINE copy_dbcsr_to_fm_sym(matrix, fm, plan)
317 : TYPE(dbcsr_type), INTENT(IN) :: matrix
318 : TYPE(cp_fm_type), INTENT(INOUT) :: fm
319 : TYPE(dbcsr_fm_copy_plan_type), INTENT(INOUT) :: plan
320 :
321 : INTEGER :: col, col_offset, i, iblock, j, k, &
322 : mirror, ncol, nrecv, nrow, pos, row, &
323 : row_offset, thread
324 : LOGICAL :: block_matches, transposed
325 : REAL(KIND=dp) :: symmetry_sign, value
326 365401 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: block
327 : TYPE(dbcsr_iterator_type) :: iter
328 : TYPE(mp_para_env_type), POINTER :: para_env
329 :
330 365353 : IF (.NOT. ALLOCATED(plan%thread)) CALL dbcsr_to_fm_plan_create(matrix, fm, plan)
331 365401 : CPASSERT(ASSOCIATED(plan%matrix_struct, fm%matrix_struct))
332 365401 : CPASSERT(plan%matrix_type == dbcsr_get_matrix_type(matrix))
333 365401 : para_env => plan%matrix_struct%context%para_env
334 365401 : symmetry_sign = 1.0_dp
335 365401 : IF (plan%matrix_type == dbcsr_type_antisymmetric) symmetry_sign = -1.0_dp
336 :
337 : !$OMP PARALLEL DEFAULT(NONE) SHARED(matrix, plan, symmetry_sign) &
338 : !$OMP PRIVATE(iter, row, col, block, transposed, row_offset, col_offset, nrow, ncol, &
339 365401 : !$OMP thread, iblock, pos, mirror, i, j, k, value, block_matches)
340 : thread = 0
341 : !$ thread = omp_get_thread_num()
342 : !$ CPASSERT(omp_get_num_threads() == SIZE(plan%thread))
343 : iblock = 0
344 : pos = 0
345 : CALL dbcsr_iterator_readonly_start(iter, matrix, dynamic=.FALSE.)
346 : DO WHILE (dbcsr_iterator_blocks_left(iter))
347 : CALL dbcsr_iterator_next_block(iter, row, col, block, transposed=transposed, &
348 : row_size=nrow, col_size=ncol, row_offset=row_offset, col_offset=col_offset)
349 : iblock = iblock + 1
350 : CPASSERT(iblock <= SIZE(plan%thread(thread)%blocks, 2))
351 : block_matches = ALL(plan%thread(thread)%blocks(:, iblock) == &
352 : [row, col, nrow, ncol, row_offset, col_offset])
353 : CPASSERT(block_matches)
354 : DO j = 1, ncol
355 : DO i = 1, nrow
356 : IF (transposed) THEN
357 : value = symmetry_sign*block(j, i)
358 : ELSE
359 : value = block(i, j)
360 : END IF
361 : DO mirror = 1, 2
362 : IF (mirror == 2 .AND. row == col) CYCLE
363 : pos = pos + 1
364 : k = plan%thread(thread)%slots(pos)
365 : plan%send_values(k) = value
366 : value = symmetry_sign*value
367 : END DO
368 : END DO
369 : END DO
370 : END DO
371 : CALL dbcsr_iterator_stop(iter)
372 : CPASSERT(iblock == SIZE(plan%thread(thread)%blocks, 2))
373 : CPASSERT(pos == SIZE(plan%thread(thread)%slots))
374 : !$OMP END PARALLEL
375 :
376 : ! Counts and destination indices were exchanged once when the plan was built.
377 : CALL para_env%alltoall(plan%send_values, plan%send_count, plan%send_disp, &
378 365401 : plan%recv_values, plan%recv_count, plan%recv_disp)
379 1093794 : nrecv = SUM(plan%recv_count)
380 193519836 : fm%local_data = 0.0_dp
381 365401 : !$OMP PARALLEL DO DEFAULT(NONE) SHARED(fm, plan, nrecv) PRIVATE(k)
382 : DO k = 1, nrecv
383 : fm%local_data(plan%recv_index(1, k), plan%recv_index(2, k)) = plan%recv_values(k)
384 : END DO
385 : !$OMP END PARALLEL DO
386 :
387 365401 : END SUBROUTINE copy_dbcsr_to_fm_sym
388 :
389 : ! **************************************************************************************************
390 : !> \brief Cache source-to-packet slots and destination FM indices, not source pointers.
391 : !> Counting and indexing use the same static DBCSR thread assignment. The
392 : !> plan is valid only while that assignment and both matrix layouts are fixed.
393 : !> \param matrix symmetric or antisymmetric DBCSR workspace
394 : !> \param fm destination full matrix
395 : !> \param plan new routing plan
396 : ! **************************************************************************************************
397 365353 : SUBROUTINE dbcsr_to_fm_plan_create(matrix, fm, plan)
398 : TYPE(dbcsr_type), INTENT(IN) :: matrix
399 : TYPE(cp_fm_type), INTENT(IN) :: fm
400 : TYPE(dbcsr_fm_copy_plan_type), INTENT(OUT) :: plan
401 :
402 : CHARACTER(LEN=*), PARAMETER :: routineN = 'dbcsr_to_fm_plan_create'
403 :
404 : INTEGER :: col, col_offset, dest, handle, i, &
405 : iblock, j, k, mirror, ncol, nproc, &
406 : nrow, nthreads, pass, pos, row, &
407 : row_offset, thread, x, y
408 : INTEGER(KIND=int_8) :: count, offset
409 365353 : INTEGER(KIND=int_8), ALLOCATABLE, DIMENSION(:, :) :: thread_pos
410 365353 : INTEGER, ALLOCATABLE, DIMENSION(:) :: nblocks
411 365353 : INTEGER, ALLOCATABLE, DIMENSION(:, :) :: send_index
412 : TYPE(cp_fm_struct_type), POINTER :: fms
413 : TYPE(dbcsr_iterator_type) :: iter
414 : TYPE(mp_para_env_type), POINTER :: para_env
415 :
416 365353 : CALL timeset(routineN, handle)
417 365353 : fms => fm%matrix_struct
418 365353 : plan%matrix_struct => fms
419 365353 : plan%matrix_type = dbcsr_get_matrix_type(matrix)
420 365353 : para_env => fms%context%para_env
421 365353 : nproc = para_env%num_pe
422 365353 : nthreads = 1
423 365353 : !$ nthreads = omp_get_max_threads()
424 2192118 : ALLOCATE (thread_pos(0:nproc - 1, 0:nthreads - 1), nblocks(0:nthreads - 1))
425 0 : ALLOCATE (plan%send_count(0:nproc - 1), plan%send_disp(0:nproc - 1), &
426 2192118 : plan%recv_count(0:nproc - 1), plan%recv_disp(0:nproc - 1))
427 365353 : thread_pos(:, :) = 0
428 365353 : nblocks(:) = 0
429 :
430 : !$OMP PARALLEL DEFAULT(NONE) &
431 : !$OMP SHARED(matrix, fms, nproc, nthreads, &
432 : !$OMP thread_pos, nblocks, send_index, plan) &
433 : !$OMP PRIVATE(iter, row, col, row_offset, col_offset, nrow, ncol, &
434 365353 : !$OMP thread, iblock, pos, pass, mirror, i, j, x, y, dest, k, offset, count)
435 : thread = 0
436 : !$ thread = omp_get_thread_num()
437 : DO pass = 1, 2
438 : iblock = 0
439 : pos = 0
440 : IF (pass == 2) THEN
441 : !$OMP SINGLE
442 : CPASSERT(SUM(thread_pos) <= INT(HUGE(0), int_8)/2)
443 : plan%send_count(:) = INT(SUM(thread_pos, DIM=2))
444 : !$ nthreads = omp_get_num_threads()
445 : ALLOCATE (plan%thread(0:nthreads - 1))
446 : DO k = 0, nthreads - 1
447 : ALLOCATE (plan%thread(k)%blocks(6, nblocks(k)), &
448 : plan%thread(k)%slots(INT(SUM(thread_pos(:, k)))))
449 : END DO
450 : offset = 0
451 : DO dest = 0, nproc - 1
452 : plan%send_disp(dest) = INT(offset)
453 : DO k = 0, nthreads - 1
454 : count = thread_pos(dest, k)
455 : thread_pos(dest, k) = offset
456 : offset = offset + count
457 : END DO
458 : END DO
459 : ALLOCATE (send_index(2, MAX(1, SUM(plan%send_count))), &
460 : plan%send_values(MAX(1, SUM(plan%send_count))))
461 : send_index(:, :) = 0
462 : plan%send_values(:) = 0.0_dp
463 : !$OMP END SINGLE
464 : END IF
465 : CALL dbcsr_iterator_readonly_start(iter, matrix, dynamic=.FALSE.)
466 : DO WHILE (dbcsr_iterator_blocks_left(iter))
467 : CALL dbcsr_iterator_next_block(iter, row, col, row_size=nrow, col_size=ncol, &
468 : row_offset=row_offset, col_offset=col_offset)
469 : iblock = iblock + 1
470 : IF (pass == 1) THEN
471 : nblocks(thread) = iblock
472 : ELSE
473 : plan%thread(thread)%blocks(:, iblock) = [row, col, nrow, ncol, row_offset, col_offset]
474 : END IF
475 : DO j = 1, ncol
476 : DO i = 1, nrow
477 : DO mirror = 1, 2
478 : IF (mirror == 2 .AND. row == col) CYCLE
479 : x = row_offset + i - 1
480 : y = col_offset + j - 1
481 : IF (mirror == 2) THEN
482 : x = col_offset + j - 1
483 : y = row_offset + i - 1
484 : END IF
485 : dest = fms%context%blacs2mpi(fms%g2p_row(x), fms%g2p_col(y))
486 : thread_pos(dest, thread) = thread_pos(dest, thread) + 1
487 : IF (pass == 1) CYCLE
488 : k = INT(thread_pos(dest, thread))
489 : send_index(:, k) = [fms%g2l_row(x), fms%g2l_col(y)]
490 : pos = pos + 1
491 : plan%thread(thread)%slots(pos) = k
492 : END DO
493 : END DO
494 : END DO
495 : END DO
496 : CALL dbcsr_iterator_stop(iter)
497 : !$OMP BARRIER
498 : END DO
499 : !$OMP END PARALLEL
500 365353 : DEALLOCATE (thread_pos, nblocks)
501 :
502 365353 : CALL para_env%alltoall(plan%send_count, plan%recv_count, 1)
503 1093650 : CPASSERT(SUM(INT(plan%recv_count, int_8)) <= INT(HUGE(0), int_8)/2)
504 365353 : plan%recv_disp(0) = 0
505 728297 : DO dest = 1, nproc - 1
506 728297 : plan%recv_disp(dest) = plan%recv_disp(dest - 1) + plan%recv_count(dest - 1)
507 : END DO
508 0 : ALLOCATE (plan%recv_index(2, MAX(1, SUM(plan%recv_count))), &
509 3283359 : plan%recv_values(MAX(1, SUM(plan%recv_count))))
510 530287930 : plan%recv_index(:, :) = 0
511 177006212 : plan%recv_values(:) = 0.0_dp
512 : CALL para_env%alltoall(send_index, 2*plan%send_count, 2*plan%send_disp, &
513 3278541 : plan%recv_index, 2*plan%recv_count, 2*plan%recv_disp)
514 365353 : CALL timestop(handle)
515 :
516 1096059 : END SUBROUTINE dbcsr_to_fm_plan_create
517 :
518 : ! **************************************************************************************************
519 : !> \brief Copy a DBCSR_BLACS matrix to a BLACS matrix
520 : !> \param bc_mat DBCSR matrix
521 : !> \param[out] fm full matrix
522 : ! **************************************************************************************************
523 887000 : SUBROUTINE copy_dbcsr_to_fm_bc(bc_mat, fm)
524 : TYPE(dbcsr_type), INTENT(IN) :: bc_mat
525 : TYPE(cp_fm_type), INTENT(INOUT) :: fm
526 :
527 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_dbcsr_to_fm_bc'
528 :
529 : INTEGER :: col, handle, row
530 887000 : INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
531 887000 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: dbcsr_block, fm_block
532 : TYPE(dbcsr_iterator_type) :: iter
533 :
534 887000 : CALL timeset(routineN, handle)
535 :
536 887000 : CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
537 :
538 : ! Now copy data to the FM matrix
539 887000 : fm_block => fm%local_data
540 202782088 : fm_block = REAL(0.0, KIND=dp)
541 : !$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
542 887000 : !$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
543 : CALL dbcsr_iterator_readonly_start(iter, bc_mat)
544 : DO WHILE (dbcsr_iterator_blocks_left(iter))
545 : CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
546 : fm_block(first_row(row):last_row(row), first_col(col):last_col(col)) = dbcsr_block(:, :)
547 : END DO
548 : CALL dbcsr_iterator_stop(iter)
549 : !$OMP END PARALLEL
550 :
551 887000 : CALL timestop(handle)
552 1774000 : END SUBROUTINE copy_dbcsr_to_fm_bc
553 :
554 : ! **************************************************************************************************
555 : !> \brief Helper routine used to copy blocks from DBCSR into FM matrices and vice versa
556 : !> \param bc_mat ...
557 : !> \param first_row ...
558 : !> \param last_row ...
559 : !> \param first_col ...
560 : !> \param last_col ...
561 : !> \author Ole Schuett
562 : ! **************************************************************************************************
563 2518060 : SUBROUTINE calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
564 : TYPE(dbcsr_type), INTENT(IN) :: bc_mat
565 : INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(OUT) :: first_row, last_row, first_col, last_col
566 :
567 : INTEGER :: col, nblkcols_local, nblkcols_total, &
568 : nblkrows_local, nblkrows_total, row
569 : INTEGER, ALLOCATABLE, DIMENSION(:) :: local_col_sizes, local_row_sizes
570 2518060 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, local_cols, local_rows, &
571 2518060 : row_blk_size
572 :
573 : CALL dbcsr_get_info(bc_mat, &
574 : nblkrows_total=nblkrows_total, &
575 : nblkcols_total=nblkcols_total, &
576 : nblkrows_local=nblkrows_local, &
577 : nblkcols_local=nblkcols_local, &
578 : local_rows=local_rows, &
579 : local_cols=local_cols, &
580 : row_blk_size=row_blk_size, &
581 2518060 : col_blk_size=col_blk_size)
582 :
583 : ! calculate first_row and last_row
584 7552924 : ALLOCATE (local_row_sizes(nblkrows_total))
585 2518060 : local_row_sizes(:) = 0
586 2518060 : IF (nblkrows_local >= 1) THEN
587 6297160 : DO row = 1, nblkrows_local
588 6297160 : local_row_sizes(local_rows(row)) = row_blk_size(local_rows(row))
589 : END DO
590 : END IF
591 7551668 : ALLOCATE (first_row(nblkrows_total), last_row(nblkrows_total))
592 2518060 : CALL dbcsr_convert_sizes_to_offsets(local_row_sizes, first_row, last_row)
593 2518060 : DEALLOCATE (local_row_sizes)
594 :
595 : ! calculate first_col and last_col
596 7550612 : ALLOCATE (local_col_sizes(nblkcols_total))
597 2518060 : local_col_sizes(:) = 0
598 2518060 : IF (nblkcols_local >= 1) THEN
599 7509411 : DO col = 1, nblkcols_local
600 7509411 : local_col_sizes(local_cols(col)) = col_blk_size(local_cols(col))
601 : END DO
602 : END IF
603 7547044 : ALLOCATE (first_col(nblkcols_total), last_col(nblkcols_total))
604 2518060 : CALL dbcsr_convert_sizes_to_offsets(local_col_sizes, first_col, last_col)
605 2518060 : DEALLOCATE (local_col_sizes)
606 :
607 2518060 : END SUBROUTINE calculate_fm_block_ranges
608 :
609 : ! **************************************************************************************************
610 : !> \brief hack for dbcsr_copy_columns
611 : !> \param matrix_b ...
612 : !> \param matrix_a ...
613 : !> \param ncol ...
614 : !> \param source_start ...
615 : !> \param target_start ...
616 : !> \param para_env ...
617 : !> \param blacs_env ...
618 : !> \author vw
619 : ! **************************************************************************************************
620 9416 : SUBROUTINE dbcsr_copy_columns_hack(matrix_b, matrix_a, &
621 : ncol, source_start, target_start, para_env, blacs_env)
622 :
623 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix_b
624 : TYPE(dbcsr_type), INTENT(IN) :: matrix_a
625 : INTEGER, INTENT(IN) :: ncol, source_start, target_start
626 : TYPE(mp_para_env_type), POINTER :: para_env
627 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
628 :
629 : INTEGER :: nfullcols_total, nfullrows_total
630 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
631 : TYPE(cp_fm_type) :: fm_matrix_a, fm_matrix_b
632 :
633 2354 : NULLIFY (fm_struct)
634 2354 : CALL dbcsr_get_info(matrix_a, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
635 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
636 2354 : ncol_global=nfullcols_total, para_env=para_env)
637 2354 : CALL cp_fm_create(fm_matrix_a, fm_struct, name="fm_matrix_a")
638 2354 : CALL cp_fm_struct_release(fm_struct)
639 :
640 2354 : CALL dbcsr_get_info(matrix_b, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
641 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
642 2354 : ncol_global=nfullcols_total, para_env=para_env)
643 2354 : CALL cp_fm_create(fm_matrix_b, fm_struct, name="fm_matrix_b")
644 2354 : CALL cp_fm_struct_release(fm_struct)
645 :
646 2354 : CALL copy_dbcsr_to_fm(matrix_a, fm_matrix_a)
647 2354 : CALL copy_dbcsr_to_fm(matrix_b, fm_matrix_b)
648 :
649 2354 : CALL cp_fm_to_fm(fm_matrix_a, fm_matrix_b, ncol, source_start, target_start)
650 :
651 2354 : CALL copy_fm_to_dbcsr(fm_matrix_b, matrix_b)
652 :
653 2354 : CALL cp_fm_release(fm_matrix_a)
654 2354 : CALL cp_fm_release(fm_matrix_b)
655 :
656 2354 : END SUBROUTINE dbcsr_copy_columns_hack
657 :
658 : ! **************************************************************************************************
659 : !> \brief Creates a DBCSR distribution from a distribution_2d
660 : !> \param[in] dist2d distribution_2d
661 : !> \param[out] dist DBCSR distribution
662 : !> \par History
663 : !> move form dbcsr_operation 01.2010
664 : ! **************************************************************************************************
665 12006 : SUBROUTINE cp_dbcsr_dist2d_to_dist(dist2d, dist)
666 : TYPE(distribution_2d_type), INTENT(IN), TARGET :: dist2d
667 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
668 :
669 12006 : INTEGER, DIMENSION(:), POINTER :: col_dist, row_dist
670 12006 : INTEGER, DIMENSION(:, :), POINTER :: col_dist_data, pgrid, row_dist_data
671 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
672 : TYPE(distribution_2d_type), POINTER :: dist2d_p
673 : TYPE(mp_para_env_type), POINTER :: para_env
674 :
675 12006 : dist2d_p => dist2d
676 : CALL distribution_2d_get(dist2d_p, &
677 : row_distribution=row_dist_data, &
678 : col_distribution=col_dist_data, &
679 12006 : blacs_env=blacs_env)
680 12006 : CALL blacs_env%get(para_env=para_env, blacs2mpi=pgrid)
681 :
682 : ! map to 1D arrays
683 12006 : row_dist => row_dist_data(:, 1)
684 12006 : col_dist => col_dist_data(:, 1)
685 : !row_cluster => row_dist_data(:, 2)
686 : !col_cluster => col_dist_data(:, 2)
687 :
688 : CALL dbcsr_distribution_new(dist, &
689 : group=para_env%get_handle(), pgrid=pgrid, &
690 : row_dist=row_dist, &
691 12006 : col_dist=col_dist)
692 :
693 12006 : END SUBROUTINE cp_dbcsr_dist2d_to_dist
694 :
695 : ! **************************************************************************************************
696 : !> \brief multiply a dbcsr with a replicated array
697 : !> c = alpha_scalar * A (dbscr) * b + c
698 : !> \param[in] matrix_a DBSCR matrxx
699 : !> \param[in] vec_b vectors b
700 : !> \param[inout] vec_c vectors c
701 : !> \param[in] ncol nbr of columns
702 : !> \param[in] alpha alpha
703 : !>
704 : ! **************************************************************************************************
705 0 : SUBROUTINE dbcsr_multiply_local(matrix_a, vec_b, vec_c, ncol, alpha)
706 : TYPE(dbcsr_type), INTENT(IN) :: matrix_a
707 : REAL(dp), DIMENSION(:, :), INTENT(IN) :: vec_b
708 : REAL(dp), DIMENSION(:, :), INTENT(INOUT) :: vec_c
709 : INTEGER, INTENT(in), OPTIONAL :: ncol
710 : REAL(dp), INTENT(IN), OPTIONAL :: alpha
711 :
712 : CHARACTER(LEN=*), PARAMETER :: routineN = 'dbcsr_multiply_local'
713 :
714 : INTEGER :: col, coloff, my_ncol, row, rowoff, &
715 : timing_handle
716 : LOGICAL :: has_symm
717 : REAL(dp) :: my_alpha, my_alpha2
718 0 : REAL(dp), DIMENSION(:, :), POINTER :: data_d
719 : TYPE(dbcsr_iterator_type) :: iter
720 :
721 0 : CALL timeset(routineN, timing_handle)
722 :
723 0 : my_alpha = 1.0_dp
724 0 : IF (PRESENT(alpha)) my_alpha = alpha
725 :
726 0 : my_ncol = SIZE(vec_b, 2)
727 0 : IF (PRESENT(ncol)) my_ncol = ncol
728 :
729 0 : my_alpha2 = 0.0_dp
730 0 : IF (dbcsr_get_matrix_type(matrix_a) == dbcsr_type_symmetric) my_alpha2 = my_alpha
731 0 : IF (dbcsr_get_matrix_type(matrix_a) == dbcsr_type_antisymmetric) my_alpha2 = -my_alpha
732 :
733 : has_symm = (dbcsr_get_matrix_type(matrix_a) == dbcsr_type_symmetric .OR. &
734 0 : dbcsr_get_matrix_type(matrix_a) == dbcsr_type_antisymmetric)
735 :
736 : !$OMP PARALLEL DEFAULT(NONE) SHARED(matrix_a,vec_b,vec_c,ncol,my_alpha2,my_alpha,my_ncol,has_symm) &
737 0 : !$OMP PRIVATE(iter,row,col,data_d,rowoff,coloff)
738 : CALL dbcsr_iterator_readonly_start(iter, matrix_a, dynamic=.TRUE., dynamic_byrows=.TRUE.)
739 : DO WHILE (dbcsr_iterator_blocks_left(iter))
740 : CALL dbcsr_iterator_next_block(iter, row, col, data_d, row_offset=rowoff, col_offset=coloff)
741 : IF (my_ncol /= 1) THEN
742 : CALL dgemm('N', 'N', &
743 : SIZE(data_d, 1), my_ncol, SIZE(data_d, 2), &
744 : my_alpha, data_d(1, 1), SIZE(data_d, 1), &
745 : vec_b(coloff, 1), SIZE(vec_b, 1), &
746 : 1.0_dp, vec_c(rowoff, 1), SIZE(vec_c, 1))
747 : ELSE
748 : CALL dgemv('N', SIZE(data_d, 1), SIZE(data_d, 2), &
749 : my_alpha, data_d(1, 1), SIZE(data_d, 1), &
750 : vec_b(coloff, 1), 1, &
751 : 1.0_dp, vec_c(rowoff, 1), 1)
752 : END IF
753 : END DO
754 : CALL dbcsr_iterator_stop(iter)
755 : !$OMP END PARALLEL
756 :
757 : ! FIXME ... in the symmetric case, the writes to vec_c depend on the column, not the row. This makes OMP-ing more difficult
758 : ! needs e.g. a buffer for vec_c and a reduction of that buffer.
759 0 : IF (has_symm) THEN
760 0 : CALL dbcsr_iterator_readonly_start(iter, matrix_a)
761 0 : DO WHILE (dbcsr_iterator_blocks_left(iter))
762 0 : CALL dbcsr_iterator_next_block(iter, row, col, data_d, row_offset=rowoff, col_offset=coloff)
763 0 : IF (row /= col) THEN
764 0 : IF (my_ncol /= 1) THEN
765 : CALL dgemm('T', 'N', &
766 : SIZE(data_d, 2), my_ncol, SIZE(data_d, 1), &
767 : my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
768 : vec_b(rowoff, 1), SIZE(vec_b, 1), &
769 0 : 1.0_dp, vec_c(coloff, 1), SIZE(vec_c, 1))
770 : ELSE
771 : CALL dgemv('T', SIZE(data_d, 1), SIZE(data_d, 2), &
772 : my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
773 : vec_b(rowoff, 1), 1, &
774 0 : 1.0_dp, vec_c(coloff, 1), 1)
775 : END IF
776 : END IF
777 : END DO
778 0 : CALL dbcsr_iterator_stop(iter)
779 : END IF
780 :
781 0 : CALL timestop(timing_handle)
782 0 : END SUBROUTINE dbcsr_multiply_local
783 :
784 : ! **************************************************************************************************
785 : !> \brief multiply a dbcsr with a fm matrix
786 : !>
787 : !> For backwards compatibility with BLAS XGEMM, this routine supports
788 : !> the multiplication of matrices with incompatible dimensions.
789 : !>
790 : !> \param[in] matrix DBCSR matrix
791 : !> \param fm_in full matrix
792 : !> \param fm_out full matrix
793 : !> \param[in] ncol nbr of columns
794 : !> \param[in] alpha alpha
795 : !> \param[in] beta beta
796 : !>
797 : ! **************************************************************************************************
798 3334278 : SUBROUTINE cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
799 : TYPE(dbcsr_type), INTENT(IN) :: matrix
800 : TYPE(cp_fm_type), INTENT(IN) :: fm_in
801 : TYPE(cp_fm_type), INTENT(INOUT) :: fm_out
802 : INTEGER, INTENT(IN) :: ncol
803 : REAL(dp), INTENT(IN), OPTIONAL :: alpha, beta
804 :
805 : CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_dbcsr_sm_fm_multiply'
806 :
807 : INTEGER :: a_ncol, a_nrow, b_ncol, b_nrow, c_ncol, &
808 : c_nrow, k_in, k_out, timing_handle, &
809 : timing_handle_mult
810 555713 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_blk_size_right_in, &
811 555713 : col_blk_size_right_out, col_dist, &
812 555713 : row_blk_size, row_dist
813 : TYPE(dbcsr_type) :: in, out
814 : TYPE(dbcsr_distribution_type) :: dist, dist_right_in, product_dist
815 : REAL(dp) :: my_alpha, my_beta
816 :
817 555713 : CALL timeset(routineN, timing_handle)
818 :
819 555713 : my_alpha = 1.0_dp
820 555713 : my_beta = 0.0_dp
821 555713 : IF (PRESENT(alpha)) my_alpha = alpha
822 555713 : IF (PRESENT(beta)) my_beta = beta
823 :
824 : ! TODO
825 555713 : CALL cp_fm_get_info(fm_in, ncol_global=b_ncol, nrow_global=b_nrow)
826 555713 : CALL cp_fm_get_info(fm_out, ncol_global=c_ncol, nrow_global=c_nrow)
827 555713 : CALL dbcsr_get_info(matrix, nfullrows_total=a_nrow, nfullcols_total=a_ncol)
828 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: A ", a_nrow, "x", a_ncol
829 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: B ", b_nrow, "x", b_ncol
830 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: C ", c_nrow, "x", c_ncol
831 :
832 555713 : CALL cp_fm_get_info(fm_out, ncol_global=k_out)
833 :
834 555713 : CALL cp_fm_get_info(fm_in, ncol_global=k_in)
835 : !write(*,*)routineN//" -----------------------------------"
836 : !IF (k_in /= k_out) &
837 : ! WRITE(*,'(3(A,I5,1X),2(A,F5.2,1X))')&
838 : ! routineN//" ncol", ncol,'k_in',k_in,'k_out',k_out,&
839 : ! 'alpha',my_alpha,'beta',my_beta
840 :
841 555713 : IF (ncol > 0 .AND. k_out > 0 .AND. k_in > 0) THEN
842 554397 : CALL dbcsr_get_info(matrix, row_blk_size=row_blk_size, col_blk_size=col_blk_size, distribution=dist)
843 554397 : CALL dbcsr_create_dist_r_unrot(dist_right_in, dist, k_in, col_blk_size_right_in)
844 :
845 : CALL dbcsr_create(in, "D", dist_right_in, dbcsr_type_no_symmetry, &
846 554397 : col_blk_size, col_blk_size_right_in)
847 :
848 554397 : CALL dbcsr_distribution_get(dist, row_dist=row_dist)
849 554397 : CALL dbcsr_distribution_get(dist_right_in, col_dist=col_dist)
850 : CALL dbcsr_distribution_new(product_dist, template=dist, &
851 554397 : row_dist=row_dist, col_dist=col_dist)
852 1663191 : ALLOCATE (col_blk_size_right_out(SIZE(col_blk_size_right_in)))
853 2254476 : col_blk_size_right_out = col_blk_size_right_in
854 554397 : CALL match_col_sizes(col_blk_size_right_out, col_blk_size_right_in, k_out)
855 :
856 : !if (k_in .ne. k_out) then
857 : ! write(*,*)routineN//" in cs", col_blk_size_right_in
858 : ! write(*,*)routineN//" out cs", col_blk_size_right_out
859 : !endif
860 :
861 : CALL dbcsr_create(out, "D", product_dist, dbcsr_type_no_symmetry, &
862 554397 : row_blk_size, col_blk_size_right_out)
863 :
864 554397 : CALL copy_fm_to_dbcsr(fm_in, in)
865 554397 : IF (ncol /= k_out .OR. my_beta /= 0.0_dp) THEN
866 126438 : CALL copy_fm_to_dbcsr(fm_out, out)
867 : END IF
868 :
869 554397 : CALL timeset(routineN//'_core', timing_handle_mult)
870 : CALL dbcsr_multiply("N", "N", my_alpha, matrix, in, my_beta, out, &
871 554397 : last_column=ncol)
872 554397 : CALL timestop(timing_handle_mult)
873 :
874 554397 : CALL copy_dbcsr_to_fm(out, fm_out)
875 :
876 554397 : CALL dbcsr_release(in)
877 554397 : CALL dbcsr_release(out)
878 554397 : DEALLOCATE (col_blk_size_right_in, col_blk_size_right_out)
879 554397 : CALL dbcsr_distribution_release(dist_right_in)
880 2217588 : CALL dbcsr_distribution_release(product_dist)
881 :
882 : END IF
883 :
884 555713 : CALL timestop(timing_handle)
885 :
886 555713 : END SUBROUTINE cp_dbcsr_sm_fm_multiply
887 :
888 : ! **************************************************************************************************
889 : !> \brief ...
890 : !> \param sizes1 ...
891 : !> \param sizes2 ...
892 : !> \param full_num ...
893 : ! **************************************************************************************************
894 554397 : SUBROUTINE match_col_sizes(sizes1, sizes2, full_num)
895 : INTEGER, DIMENSION(:), INTENT(INOUT) :: sizes1
896 : INTEGER, DIMENSION(:), INTENT(IN) :: sizes2
897 : INTEGER, INTENT(IN) :: full_num
898 :
899 : INTEGER :: left, n1, n2, p, rm, used
900 :
901 554397 : n1 = SIZE(sizes1)
902 554397 : n2 = SIZE(sizes2)
903 554397 : IF (n1 /= n2) THEN
904 0 : CPABORT("distributions must be equal!")
905 : END IF
906 1127238 : sizes1(1:n1) = sizes2(1:n1)
907 1127238 : used = SUM(sizes1(1:n1))
908 : ! If sizes1 does not cover everything, then we increase the
909 : ! size of the last block; otherwise we reduce the blocks
910 : ! (from the end) until it is small enough.
911 554397 : IF (used < full_num) THEN
912 0 : sizes1(n1) = sizes1(n1) + full_num - used
913 : ELSE
914 554397 : left = used - full_num
915 554397 : p = n1
916 554397 : DO WHILE (left > 0 .AND. p > 0)
917 0 : rm = MIN(left, sizes1(p))
918 0 : sizes1(p) = sizes1(p) - rm
919 0 : left = left - rm
920 0 : p = p - 1
921 : END DO
922 : END IF
923 554397 : END SUBROUTINE match_col_sizes
924 :
925 : ! **************************************************************************************************
926 : !> \brief performs the multiplication sparse_matrix+dense_mat*dens_mat^T
927 : !> if matrix_g is not explicitly given, matrix_v^T will be used
928 : !> this can be important to save the necessary redistribute for a
929 : !> different matrix_g and increase performance.
930 : !> \param sparse_matrix ...
931 : !> \param matrix_v ...
932 : !> \param matrix_g ...
933 : !> \param ncol ...
934 : !> \param alpha ...
935 : !> \param keep_sparsity Determines if the sparsity of sparse_matrix is retained
936 : !> by default it is TRUE
937 : !> \param symmetry_mode There are the following modes
938 : !> 1: sparse_matrix += 0.5*alpha*(v*g^T+g^T*v) (symmetric update)
939 : !> -1: sparse_matrix += 0.5*alpha*(v*g^T-g^T*v) (skewsymmetric update)
940 : !> else: sparse_matrix += alpha*v*g^T (no symmetry, default)
941 : !> saves some redistribution steps
942 : ! **************************************************************************************************
943 262854 : SUBROUTINE cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
944 : TYPE(dbcsr_type), INTENT(INOUT) :: sparse_matrix
945 : TYPE(cp_fm_type), INTENT(IN) :: matrix_v
946 : TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: matrix_g
947 : INTEGER, INTENT(IN) :: ncol
948 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: alpha
949 : LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
950 : INTEGER, INTENT(IN), OPTIONAL :: symmetry_mode
951 :
952 : CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_dbcsr_plus_fm_fm_t'
953 :
954 : INTEGER :: k, my_symmetry_mode, nao, npcols, &
955 : timing_handle
956 262854 : INTEGER, DIMENSION(:), POINTER :: col_blk_size_left, col_dist_left, &
957 262854 : row_blk_size, row_dist
958 : LOGICAL :: check_product, my_keep_sparsity
959 : REAL(KIND=dp) :: my_alpha, norm
960 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
961 : TYPE(cp_fm_type) :: fm_matrix
962 : TYPE(dbcsr_distribution_type) :: dist_left, sparse_dist
963 : TYPE(dbcsr_type) :: mat_g, mat_v, sparse_matrix2, &
964 : sparse_matrix3
965 :
966 262854 : check_product = .FALSE.
967 :
968 262854 : CALL timeset(routineN, timing_handle)
969 :
970 262854 : my_keep_sparsity = .TRUE.
971 262854 : IF (PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
972 :
973 262854 : my_symmetry_mode = 0
974 262854 : IF (PRESENT(symmetry_mode)) my_symmetry_mode = symmetry_mode
975 :
976 262854 : NULLIFY (col_dist_left)
977 :
978 262854 : IF (ncol > 0) THEN
979 260778 : IF (.NOT. dbcsr_valid_index(sparse_matrix)) THEN
980 0 : CPABORT("sparse_matrix must pre-exist")
981 : END IF
982 : !
983 : ! Setup matrix_v
984 260778 : CALL cp_fm_get_info(matrix_v, ncol_global=k)
985 : !WRITE(*,*)routineN//'truncated mult k, ncol',k,ncol,' PRESENT (matrix_g)',PRESENT (matrix_g)
986 260778 : CALL dbcsr_get_info(sparse_matrix, distribution=sparse_dist)
987 260778 : CALL dbcsr_distribution_get(sparse_dist, npcols=npcols, row_dist=row_dist)
988 260778 : CALL create_bl_distribution(col_dist_left, col_blk_size_left, k, npcols)
989 : CALL dbcsr_distribution_new(dist_left, template=sparse_dist, &
990 260778 : row_dist=row_dist, col_dist=col_dist_left)
991 260778 : DEALLOCATE (col_dist_left)
992 260778 : CALL dbcsr_get_info(sparse_matrix, row_blk_size=row_blk_size)
993 : CALL dbcsr_create(mat_v, "DBCSR matrix_v", dist_left, dbcsr_type_no_symmetry, &
994 260778 : row_blk_size, col_blk_size_left)
995 260778 : CALL copy_fm_to_dbcsr(matrix_v, mat_v)
996 260778 : CALL dbcsr_verify_matrix(mat_v)
997 : !
998 : ! Setup matrix_g
999 260778 : IF (PRESENT(matrix_g)) THEN
1000 : CALL dbcsr_create(mat_g, "DBCSR matrix_g", dist_left, dbcsr_type_no_symmetry, &
1001 127397 : row_blk_size, col_blk_size_left)
1002 127397 : CALL copy_fm_to_dbcsr(matrix_g, mat_g)
1003 : END IF
1004 : !
1005 260778 : DEALLOCATE (col_blk_size_left)
1006 260778 : CALL dbcsr_distribution_release(dist_left)
1007 : !
1008 : !
1009 : IF (check_product) THEN
1010 : CALL cp_fm_get_info(matrix_v, nrow_global=nao)
1011 : CALL cp_fm_struct_create(fm_struct_tmp, context=matrix_v%matrix_struct%context, nrow_global=nao, &
1012 : ncol_global=nao, para_env=matrix_v%matrix_struct%para_env)
1013 : CALL cp_fm_create(fm_matrix, fm_struct_tmp, name="fm matrix")
1014 : CALL cp_fm_struct_release(fm_struct_tmp)
1015 : CALL copy_dbcsr_to_fm(sparse_matrix, fm_matrix)
1016 : CALL dbcsr_copy(sparse_matrix3, sparse_matrix)
1017 : END IF
1018 : !
1019 260778 : my_alpha = 1.0_dp
1020 260778 : IF (PRESENT(alpha)) my_alpha = alpha
1021 260778 : IF (PRESENT(matrix_g)) THEN
1022 127397 : IF (my_symmetry_mode == 1) THEN
1023 : ! Symmetric mode
1024 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
1025 : 1.0_dp, sparse_matrix, &
1026 : retain_sparsity=my_keep_sparsity, &
1027 44152 : last_k=ncol)
1028 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_g, mat_v, &
1029 : 1.0_dp, sparse_matrix, &
1030 : retain_sparsity=my_keep_sparsity, &
1031 44152 : last_k=ncol)
1032 83245 : ELSE IF (my_symmetry_mode == -1) THEN
1033 : ! Skewsymmetric mode
1034 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
1035 : 1.0_dp, sparse_matrix, &
1036 : retain_sparsity=my_keep_sparsity, &
1037 2594 : last_k=ncol)
1038 : CALL dbcsr_multiply("N", "T", -0.5_dp*my_alpha, mat_g, mat_v, &
1039 : 1.0_dp, sparse_matrix, &
1040 : retain_sparsity=my_keep_sparsity, &
1041 2594 : last_k=ncol)
1042 : ELSE
1043 : ! Normal mode
1044 : CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_g, &
1045 : 1.0_dp, sparse_matrix, &
1046 : retain_sparsity=my_keep_sparsity, &
1047 80651 : last_k=ncol)
1048 : END IF
1049 : ELSE
1050 : CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_v, &
1051 : 1.0_dp, sparse_matrix, &
1052 : retain_sparsity=my_keep_sparsity, &
1053 133381 : last_k=ncol)
1054 : END IF
1055 :
1056 : IF (check_product) THEN
1057 : IF (PRESENT(matrix_g)) THEN
1058 : IF (my_symmetry_mode == 1) THEN
1059 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
1060 : 1.0_dp, fm_matrix)
1061 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_g, matrix_v, &
1062 : 1.0_dp, fm_matrix)
1063 : ELSE IF (my_symmetry_mode == -1) THEN
1064 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
1065 : 1.0_dp, fm_matrix)
1066 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, -0.5_dp*my_alpha, matrix_g, matrix_v, &
1067 : 1.0_dp, fm_matrix)
1068 : ELSE
1069 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_g, &
1070 : 1.0_dp, fm_matrix)
1071 : END IF
1072 : ELSE
1073 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_v, &
1074 : 1.0_dp, fm_matrix)
1075 : END IF
1076 :
1077 : CALL dbcsr_copy(sparse_matrix2, sparse_matrix)
1078 : CALL dbcsr_scale(sparse_matrix2, alpha_scalar=0.0_dp)
1079 : CALL copy_fm_to_dbcsr(fm_matrix, sparse_matrix2, keep_sparsity=my_keep_sparsity)
1080 : CALL dbcsr_add(sparse_matrix2, sparse_matrix, alpha_scalar=1.0_dp, &
1081 : beta_scalar=-1.0_dp)
1082 : norm = dbcsr_frobenius_norm(sparse_matrix2)
1083 : WRITE (*, *) 'nao=', nao, ' k=', k, ' ncol=', ncol, ' my_alpha=', my_alpha
1084 : WRITE (*, *) 'PRESENT (matrix_g)', PRESENT(matrix_g)
1085 : WRITE (*, *) 'matrix_type=', dbcsr_get_matrix_type(sparse_matrix)
1086 : WRITE (*, *) 'norm(sm+alpha*v*g^t - fm+alpha*v*g^t)/n=', norm/REAL(nao, dp)
1087 : CALL dbcsr_release(sparse_matrix2)
1088 : CALL dbcsr_release(sparse_matrix3)
1089 : CALL cp_fm_release(fm_matrix)
1090 : END IF
1091 260778 : CALL dbcsr_release(mat_v)
1092 260778 : IF (PRESENT(matrix_g)) CALL dbcsr_release(mat_g)
1093 : END IF
1094 262854 : CALL timestop(timing_handle)
1095 :
1096 262854 : END SUBROUTINE cp_dbcsr_plus_fm_fm_t
1097 :
1098 : ! **************************************************************************************************
1099 : !> \brief Utility function to copy a specially shaped fm to dbcsr_matrix
1100 : !> The result matrix will be the matrix in dbcsr format
1101 : !> with the row blocks sizes according to the block_sizes of the template
1102 : !> and the col blocks sizes evenly blocked with the internal dbcsr conversion
1103 : !> size (32 is the current default)
1104 : !> \param matrix ...
1105 : !> \param fm_in ...
1106 : !> \param template ...
1107 : ! **************************************************************************************************
1108 16785 : SUBROUTINE cp_fm_to_dbcsr_row_template(matrix, fm_in, template)
1109 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix
1110 : TYPE(cp_fm_type), INTENT(IN) :: fm_in
1111 : TYPE(dbcsr_type), INTENT(IN) :: template
1112 :
1113 : INTEGER :: k_in
1114 5595 : INTEGER, DIMENSION(:), POINTER :: col_blk_size_right_in, row_blk_size
1115 : TYPE(dbcsr_distribution_type) :: dist_right_in, tmpl_dist
1116 :
1117 5595 : CALL cp_fm_get_info(fm_in, ncol_global=k_in)
1118 :
1119 5595 : CALL dbcsr_get_info(template, distribution=tmpl_dist)
1120 5595 : CALL dbcsr_create_dist_r_unrot(dist_right_in, tmpl_dist, k_in, col_blk_size_right_in)
1121 5595 : CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
1122 : CALL dbcsr_create(matrix, "D", dist_right_in, dbcsr_type_no_symmetry, &
1123 5595 : row_blk_size, col_blk_size_right_in)
1124 :
1125 5595 : CALL copy_fm_to_dbcsr(fm_in, matrix)
1126 5595 : DEALLOCATE (col_blk_size_right_in)
1127 5595 : CALL dbcsr_distribution_release(dist_right_in)
1128 :
1129 5595 : END SUBROUTINE cp_fm_to_dbcsr_row_template
1130 :
1131 : ! **************************************************************************************************
1132 : !> \brief Utility function to create an arbitrary shaped dbcsr matrix
1133 : !> with the same processor grid as the template matrix
1134 : !> both row sizes and col sizes are evenly blocked with the internal
1135 : !> dbcsr_conversion size (32 is the current default)
1136 : !> \param matrix dbcsr matrix to be created
1137 : !> \param template template dbcsr matrix giving its mp_env
1138 : !> \param m global row size of output matrix
1139 : !> \param n global col size of output matrix
1140 : !> \param sym ...
1141 : ! **************************************************************************************************
1142 347380 : SUBROUTINE cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym)
1143 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
1144 : INTEGER, INTENT(IN) :: m, n
1145 : CHARACTER, INTENT(IN), OPTIONAL :: sym
1146 :
1147 : CHARACTER :: mysym
1148 : INTEGER :: npcols, nprows
1149 173690 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
1150 173690 : row_dist
1151 : TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
1152 :
1153 173690 : CALL dbcsr_get_info(template, matrix_type=mysym, distribution=tmpl_dist)
1154 :
1155 173690 : IF (PRESENT(sym)) mysym = sym
1156 :
1157 173690 : NULLIFY (row_dist, col_dist)
1158 173690 : NULLIFY (row_blk_size, col_blk_size)
1159 : !NULLIFY (row_cluster, col_cluster)
1160 :
1161 173690 : CALL dbcsr_distribution_get(tmpl_dist, nprows=nprows, npcols=npcols)
1162 173690 : CALL create_bl_distribution(row_dist, row_blk_size, m, nprows)
1163 173690 : CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
1164 : CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
1165 : row_dist=row_dist, col_dist=col_dist, &
1166 : !row_cluster=row_cluster, col_cluster=col_cluster, &
1167 173690 : reuse_arrays=.TRUE.)
1168 :
1169 : CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, &
1170 173690 : row_blk_size, col_blk_size, reuse_arrays=.TRUE.)
1171 173690 : CALL dbcsr_distribution_release(dist_m_n)
1172 :
1173 173690 : END SUBROUTINE cp_dbcsr_m_by_n_from_template
1174 :
1175 : ! **************************************************************************************************
1176 : !> \brief Utility function to create dbcsr matrix, m x n matrix (n arbitrary)
1177 : !> with the same processor grid and row distribution as the template matrix
1178 : !> col sizes are evenly blocked with the internal
1179 : !> dbcsr_conversion size (32 is the current default)
1180 : !> \param matrix dbcsr matrix to be created
1181 : !> \param template template dbcsr matrix giving its mp_env
1182 : !> \param n global col size of output matrix
1183 : !> \param sym ...
1184 : ! **************************************************************************************************
1185 652688 : SUBROUTINE cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym)
1186 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
1187 : INTEGER :: n
1188 : CHARACTER, OPTIONAL :: sym
1189 :
1190 : CHARACTER :: mysym
1191 : INTEGER :: npcols
1192 163172 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
1193 163172 : row_dist
1194 : TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
1195 :
1196 326344 : mysym = dbcsr_get_matrix_type(template)
1197 163172 : IF (PRESENT(sym)) mysym = sym
1198 :
1199 163172 : CALL dbcsr_get_info(template, distribution=tmpl_dist)
1200 : CALL dbcsr_distribution_get(tmpl_dist, &
1201 : npcols=npcols, &
1202 163172 : row_dist=row_dist)
1203 :
1204 163172 : NULLIFY (col_dist, col_blk_size)
1205 163172 : CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
1206 : CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
1207 163172 : row_dist=row_dist, col_dist=col_dist)
1208 :
1209 163172 : CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
1210 163172 : CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, row_blk_size, col_blk_size)
1211 :
1212 163172 : DEALLOCATE (col_dist, col_blk_size)
1213 163172 : CALL dbcsr_distribution_release(dist_m_n)
1214 :
1215 163172 : END SUBROUTINE cp_dbcsr_m_by_n_from_row_template
1216 :
1217 : ! **************************************************************************************************
1218 : !> \brief Distributes elements into blocks and into bins
1219 : !>
1220 : !> \param[out] block_distribution block distribution to bins
1221 : !> \param[out] block_size sizes of blocks
1222 : !> \param[in] nelements number of elements to bin
1223 : !> \param[in] nbins number of bins
1224 : !> \par Term clarification
1225 : !> An example: blocks are atom blocks and bins are process rows/columns.
1226 : ! **************************************************************************************************
1227 1331322 : SUBROUTINE create_bl_distribution(block_distribution, &
1228 : block_size, nelements, nbins)
1229 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: block_distribution, block_size
1230 : INTEGER, INTENT(IN) :: nelements, nbins
1231 :
1232 : CHARACTER(len=*), PARAMETER :: routineN = 'create_bl_distribution', &
1233 : routineP = moduleN//':'//routineN
1234 :
1235 : INTEGER :: bin, blk_layer, element_stack, els, &
1236 : estimated_blocks, max_blocks_per_bin, &
1237 : nblks, nblocks, stat
1238 1331322 : INTEGER, DIMENSION(:), POINTER :: blk_dist, blk_sizes
1239 :
1240 : ! ---------------------------------------------------------------------------
1241 :
1242 1331322 : NULLIFY (block_distribution)
1243 1331322 : NULLIFY (block_size)
1244 : ! Define the sizes on which we build the distribution.
1245 1331322 : IF (nelements > 0) THEN
1246 :
1247 1316560 : nblocks = CEILING(REAL(nelements, KIND=dp)/REAL(max_elements_per_block, KIND=dp))
1248 1316560 : max_blocks_per_bin = CEILING(REAL(nblocks, KIND=dp)/REAL(nbins, KIND=dp))
1249 :
1250 : IF (debug_mod) THEN
1251 : WRITE (*, '(1X,A,1X,A,I7,A,I7,A)') routineP, "For", nelements, &
1252 : " elements and", nbins, " bins"
1253 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1254 : max_elements_per_block, " max elements per block"
1255 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1256 : nblocks, " blocks"
1257 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1258 : max_blocks_per_bin, " max blocks/bin"
1259 : END IF
1260 :
1261 1316560 : estimated_blocks = max_blocks_per_bin*nbins
1262 3949680 : ALLOCATE (blk_dist(estimated_blocks), stat=stat)
1263 1316560 : IF (stat /= 0) THEN
1264 0 : CPABORT("blk_dist")
1265 : END IF
1266 2633120 : ALLOCATE (blk_sizes(estimated_blocks), stat=stat)
1267 : IF (stat /= 0) THEN
1268 0 : CPABORT("blk_sizes")
1269 : END IF
1270 1316560 : element_stack = 0
1271 1316560 : nblks = 0
1272 2711384 : DO blk_layer = 1, max_blocks_per_bin
1273 4255694 : DO bin = 0, nbins - 1
1274 1544310 : els = MIN(max_elements_per_block, nelements - element_stack)
1275 2939134 : IF (els > 0) THEN
1276 1405680 : element_stack = element_stack + els
1277 1405680 : nblks = nblks + 1
1278 1405680 : blk_dist(nblks) = bin
1279 1405680 : blk_sizes(nblks) = els
1280 : IF (debug_mod) WRITE (*, '(1X,A,I5,A,I5,A,I5)') routineP//" Assigning", &
1281 : els, " elements as block", nblks, " to bin", bin
1282 : END IF
1283 : END DO
1284 : END DO
1285 : ! Create the output arrays.
1286 1316560 : IF (nblks == estimated_blocks) THEN
1287 1177930 : block_distribution => blk_dist
1288 1177930 : block_size => blk_sizes
1289 : ELSE
1290 415890 : ALLOCATE (block_distribution(nblks), stat=stat)
1291 : IF (stat /= 0) THEN
1292 0 : CPABORT("blk_dist")
1293 : END IF
1294 557800 : block_distribution(:) = blk_dist(1:nblks)
1295 138630 : DEALLOCATE (blk_dist)
1296 277260 : ALLOCATE (block_size(nblks), stat=stat)
1297 : IF (stat /= 0) THEN
1298 0 : CPABORT("blk_sizes")
1299 : END IF
1300 557800 : block_size(:) = blk_sizes(1:nblks)
1301 138630 : DEALLOCATE (blk_sizes)
1302 : END IF
1303 : ELSE
1304 14762 : ALLOCATE (block_distribution(0), stat=stat)
1305 : IF (stat /= 0) THEN
1306 0 : CPABORT("blk_dist")
1307 : END IF
1308 14762 : ALLOCATE (block_size(0), stat=stat)
1309 : IF (stat /= 0) THEN
1310 0 : CPABORT("blk_sizes")
1311 : END IF
1312 : END IF
1313 : 1579 FORMAT(I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5)
1314 : IF (debug_mod) THEN
1315 : WRITE (*, '(1X,A,A)') routineP//" Distribution"
1316 : WRITE (*, 1579) block_distribution(:)
1317 : WRITE (*, '(1X,A,A)') routineP//" Sizes"
1318 : WRITE (*, 1579) block_size(:)
1319 : END IF
1320 1331322 : END SUBROUTINE create_bl_distribution
1321 :
1322 : ! **************************************************************************************************
1323 : !> \brief Creates a new distribution for the right matrix in a matrix
1324 : !> multiplication with unrotated grid.
1325 : !> \param[out] dist_right new distribution for the right matrix
1326 : !> \param[in] dist_left the distribution of the left matrix
1327 : !> \param[in] ncolumns number of columns in right matrix
1328 : !> \param[out] right_col_blk_sizes sizes of blocks in the created column
1329 : !> \par The new row distribution for the right matrix is the same as the row
1330 : !> distribution of the left matrix, while the column distribution is
1331 : !> created so that it is appropriate to the parallel environment.
1332 : ! **************************************************************************************************
1333 559992 : SUBROUTINE dbcsr_create_dist_r_unrot(dist_right, dist_left, ncolumns, &
1334 : right_col_blk_sizes)
1335 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist_right
1336 : TYPE(dbcsr_distribution_type), INTENT(IN) :: dist_left
1337 : INTEGER, INTENT(IN) :: ncolumns
1338 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: right_col_blk_sizes
1339 :
1340 : INTEGER :: multiplicity, ncols, nimages, npcols, &
1341 : nprows
1342 : INTEGER, ALLOCATABLE, DIMENSION(:) :: tmp_images
1343 559992 : INTEGER, DIMENSION(:), POINTER :: old_col_dist, right_col_dist, &
1344 559992 : right_row_dist
1345 :
1346 : CALL dbcsr_distribution_get(dist_left, &
1347 : ncols=ncols, &
1348 : col_dist=old_col_dist, &
1349 : nprows=nprows, &
1350 559992 : npcols=npcols)
1351 :
1352 : ! Create the column distribution
1353 559992 : CALL create_bl_distribution(right_col_dist, right_col_blk_sizes, ncolumns, npcols)
1354 : ! Create an even row distribution.
1355 2239968 : ALLOCATE (right_row_dist(ncols), tmp_images(ncols))
1356 559992 : nimages = lcm(nprows, npcols)/nprows
1357 559992 : multiplicity = nprows/gcd(nprows, npcols)
1358 559992 : CALL rebin_distribution(right_row_dist, tmp_images, old_col_dist, nprows, multiplicity, nimages)
1359 :
1360 : CALL dbcsr_distribution_new(dist_right, &
1361 : template=dist_left, &
1362 : row_dist=right_row_dist, &
1363 : col_dist=right_col_dist, &
1364 : !row_cluster=dummy,&
1365 : !col_cluster=dummy,&
1366 559992 : reuse_arrays=.TRUE.)
1367 559992 : DEALLOCATE (tmp_images)
1368 559992 : END SUBROUTINE dbcsr_create_dist_r_unrot
1369 :
1370 : ! **************************************************************************************************
1371 : !> \brief Makes new distribution with decimation and multiplicity
1372 : !> \param[out] new_bins new real distribution
1373 : !> \param[out] images new image distribution
1374 : !> \param[in] source_bins Basis for the new distribution and images
1375 : !> \param[in] nbins number of bins in the new real distribution
1376 : !> \param[in] multiplicity multiplicity
1377 : !> \param[in] nimages number of images in the new distribution
1378 : !> \par Definition of multiplicity and nimages
1379 : !> Multiplicity and decimation (number of images) are used to
1380 : !> match process grid coordinates on non-square process
1381 : !> grids. Given source_nbins and target_nbins, their relation is
1382 : !> source_nbins * target_multiplicity
1383 : !> = target_nbins * target_nimages.
1384 : !> It is best when both multiplicity and nimages are small. To
1385 : !> get these two factors, then, one can use the following formulas:
1386 : !> nimages = lcm(source_nbins, target_nbins) / target_nbins
1387 : !> multiplicity = target_nbins / gcd(source_nbins, target_nbins)
1388 : !> from the target's point of view (nimages = target_nimages).
1389 : !> \par Mapping
1390 : !> The new distribution comprises of real bins and images within
1391 : !> bins. These can be view as target_nbins*nimages virtual
1392 : !> columns. These same virtual columns are also
1393 : !> source_nbins*multiplicity in number. Therefore these virtual
1394 : !> columns are mapped from source_nbins*multiplicity onto
1395 : !> target_bins*nimages (each target bin has nimages images):
1396 : !> Source 4: |1 2 3|4 5 6|7 8 9|A B C| (4*3)
1397 : !> Target 6: |1 2|3 4|5 6|7 8|9 A|B C| (6*2)
1398 : !> multiplicity=3, nimages=2, 12 virtual columns (1-C).
1399 : !> Source bin elements are evenly mapped into one of multiplicity
1400 : !> virtual columns. Other (non-even, block-size aware) mappings
1401 : !> could be better.
1402 : ! **************************************************************************************************
1403 559992 : SUBROUTINE rebin_distribution(new_bins, images, source_bins, &
1404 : nbins, multiplicity, nimages)
1405 : INTEGER, DIMENSION(:), INTENT(OUT) :: new_bins, images
1406 : INTEGER, DIMENSION(:), INTENT(IN) :: source_bins
1407 : INTEGER, INTENT(IN) :: nbins, multiplicity, nimages
1408 :
1409 : INTEGER :: bin, i, old_nbins, virtual_bin
1410 559992 : INTEGER, ALLOCATABLE, DIMENSION(:) :: bin_multiplier
1411 :
1412 : ! ---------------------------------------------------------------------------
1413 :
1414 559992 : IF (MOD(nbins*nimages, multiplicity) /= 0) THEN
1415 0 : CPWARN("mulitplicity is not divisor of new process grid coordinate")
1416 : END IF
1417 559992 : old_nbins = (nbins*nimages)/multiplicity
1418 1679976 : ALLOCATE (bin_multiplier(0:old_nbins - 1))
1419 559992 : bin_multiplier(:) = 0
1420 2600643 : DO i = 1, SIZE(new_bins)
1421 2040651 : IF (i <= SIZE(source_bins)) THEN
1422 2040651 : bin = source_bins(i)
1423 : ELSE
1424 : ! Fill remainder with a cyclic distribution
1425 0 : bin = MOD(i, old_nbins)
1426 : END IF
1427 2040651 : virtual_bin = bin*multiplicity + bin_multiplier(bin)
1428 2040651 : new_bins(i) = virtual_bin/nimages
1429 2040651 : images(i) = 1 + MOD(virtual_bin, nimages)
1430 2040651 : bin_multiplier(bin) = bin_multiplier(bin) + 1
1431 2600643 : IF (bin_multiplier(bin) >= multiplicity) THEN
1432 874193 : bin_multiplier(bin) = 0
1433 : END IF
1434 : END DO
1435 559992 : END SUBROUTINE rebin_distribution
1436 :
1437 : ! **************************************************************************************************
1438 : !> \brief Creates a block-cyclic compatible distribution
1439 : !>
1440 : !> All blocks in a dimension, except for possibly the last
1441 : !> block, have the same size.
1442 : !> \param[out] dist the elemental distribution
1443 : !> \param[in] nrows number of full rows
1444 : !> \param[in] ncolumns number of full columns
1445 : !> \param[in] nrow_block size of row blocks
1446 : !> \param[in] ncol_block size of column blocks
1447 : !> \param group_handle ...
1448 : !> \param pgrid ...
1449 : !> \param[out] row_blk_sizes row block sizes
1450 : !> \param[out] col_blk_sizes column block sizes
1451 : ! **************************************************************************************************
1452 2518060 : SUBROUTINE dbcsr_create_dist_block_cyclic(dist, nrows, ncolumns, &
1453 : nrow_block, ncol_block, group_handle, pgrid, row_blk_sizes, col_blk_sizes)
1454 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
1455 : INTEGER, INTENT(IN) :: nrows, ncolumns, nrow_block, ncol_block, &
1456 : group_handle
1457 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
1458 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: row_blk_sizes, col_blk_sizes
1459 :
1460 : CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_create_dist_block_cyclic'
1461 :
1462 : INTEGER :: nblkcols, nblkrows, npcols, nprows, &
1463 : pdim, sz
1464 2518060 : INTEGER, DIMENSION(:), POINTER :: cd_data, rd_data
1465 :
1466 : ! Row sizes
1467 2518060 : IF (nrow_block == 0) THEN
1468 : nblkrows = 0
1469 : sz = 0
1470 : ELSE
1471 2518060 : nblkrows = nrows/nrow_block
1472 2518060 : sz = MOD(nrows, nrow_block)
1473 : END IF
1474 2518060 : IF (sz > 0) nblkrows = nblkrows + 1
1475 10069728 : ALLOCATE (row_blk_sizes(nblkrows), rd_data(nblkrows))
1476 9909432 : row_blk_sizes = nrow_block
1477 2518060 : IF (sz /= 0) row_blk_sizes(nblkrows) = sz
1478 :
1479 : ! Column sizes
1480 2518060 : IF (ncol_block == 0) THEN
1481 : nblkcols = 0
1482 : sz = 0
1483 : ELSE
1484 2518060 : nblkcols = ncolumns/ncol_block
1485 2518060 : sz = MOD(ncolumns, ncol_block)
1486 : END IF
1487 2518060 : IF (sz > 0) nblkcols = nblkcols + 1
1488 10065104 : ALLOCATE (col_blk_sizes(nblkcols), cd_data(nblkcols))
1489 7512979 : col_blk_sizes = ncol_block
1490 2518060 : IF (sz /= 0) col_blk_sizes(nblkcols) = sz
1491 : !
1492 : IF (debug_mod) THEN
1493 : WRITE (*, *) routineN//" nrows,nrow_block,nblkrows=", &
1494 : nrows, nrow_block, nblkrows
1495 : WRITE (*, *) routineN//" ncols,ncol_block,nblkcols=", &
1496 : ncolumns, ncol_block, nblkcols
1497 : END IF
1498 : ! Calculate process row distribution
1499 2518060 : nprows = SIZE(pgrid, 1)
1500 7420852 : DO pdim = 0, MIN(nprows - 1, nblkrows - 1)
1501 14812224 : rd_data(1 + pdim:nblkrows:nprows) = pdim
1502 : END DO
1503 : ! Calculate process column distribution
1504 2518060 : npcols = SIZE(pgrid, 2)
1505 5032552 : DO pdim = 0, MIN(npcols - 1, nblkcols - 1)
1506 10027471 : cd_data(1 + pdim:nblkcols:npcols) = pdim
1507 : END DO
1508 : !
1509 : IF (debug_mod) THEN
1510 : WRITE (*, *) routineN//" row_dist", &
1511 : rd_data
1512 : WRITE (*, *) routineN//" col_dist", &
1513 : cd_data
1514 : END IF
1515 : !
1516 : CALL dbcsr_distribution_new(dist, &
1517 : group=group_handle, pgrid=pgrid, &
1518 : row_dist=rd_data, &
1519 : col_dist=cd_data, &
1520 2518060 : reuse_arrays=.TRUE.)
1521 :
1522 2518060 : END SUBROUTINE dbcsr_create_dist_block_cyclic
1523 :
1524 : ! **************************************************************************************************
1525 : !> \brief Allocate and initialize a real matrix 1-dimensional set.
1526 : !> \param[in,out] matrix_set Set containing the DBCSR matrices
1527 : !> \param[in] nmatrix Size of set
1528 : !> \par History
1529 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1530 : ! **************************************************************************************************
1531 204495 : SUBROUTINE allocate_dbcsr_matrix_set_1d(matrix_set, nmatrix)
1532 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1533 : INTEGER, INTENT(IN) :: nmatrix
1534 :
1535 : INTEGER :: imatrix
1536 :
1537 204495 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1538 1415284 : ALLOCATE (matrix_set(nmatrix))
1539 1006294 : DO imatrix = 1, nmatrix
1540 1006294 : NULLIFY (matrix_set(imatrix)%matrix)
1541 : END DO
1542 204495 : END SUBROUTINE allocate_dbcsr_matrix_set_1d
1543 :
1544 : ! **************************************************************************************************
1545 : !> \brief Allocate and initialize a real matrix 2-dimensional set.
1546 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1547 : !> \param[in] nmatrix Size of set
1548 : !> \param mmatrix ...
1549 : !> \par History
1550 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1551 : ! **************************************************************************************************
1552 195676 : SUBROUTINE allocate_dbcsr_matrix_set_2d(matrix_set, nmatrix, mmatrix)
1553 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1554 : INTEGER, INTENT(IN) :: nmatrix, mmatrix
1555 :
1556 : INTEGER :: imatrix, jmatrix
1557 :
1558 195676 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1559 5026416 : ALLOCATE (matrix_set(nmatrix, mmatrix))
1560 2186988 : DO jmatrix = 1, mmatrix
1561 4439388 : DO imatrix = 1, nmatrix
1562 4243712 : NULLIFY (matrix_set(imatrix, jmatrix)%matrix)
1563 : END DO
1564 : END DO
1565 195676 : END SUBROUTINE allocate_dbcsr_matrix_set_2d
1566 :
1567 : ! **************************************************************************************************
1568 : !> \brief Allocate and initialize a real matrix 3-dimensional set.
1569 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1570 : !> \param[in] nmatrix Size of set
1571 : !> \param mmatrix ...
1572 : !> \param pmatrix ...
1573 : !> \par History
1574 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1575 : ! **************************************************************************************************
1576 328 : SUBROUTINE allocate_dbcsr_matrix_set_3d(matrix_set, nmatrix, mmatrix, pmatrix)
1577 : TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1578 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix
1579 :
1580 : INTEGER :: imatrix, jmatrix, kmatrix
1581 :
1582 328 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1583 7662 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix))
1584 668 : DO kmatrix = 1, pmatrix
1585 2562 : DO jmatrix = 1, mmatrix
1586 6022 : DO imatrix = 1, nmatrix
1587 5682 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1588 : END DO
1589 : END DO
1590 : END DO
1591 328 : END SUBROUTINE allocate_dbcsr_matrix_set_3d
1592 :
1593 : ! **************************************************************************************************
1594 : !> \brief Allocate and initialize a real matrix 4-dimensional set.
1595 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1596 : !> \param[in] nmatrix Size of set
1597 : !> \param mmatrix ...
1598 : !> \param pmatrix ...
1599 : !> \param qmatrix ...
1600 : !> \par History
1601 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1602 : ! **************************************************************************************************
1603 0 : SUBROUTINE allocate_dbcsr_matrix_set_4d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix)
1604 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1605 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix
1606 :
1607 : INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1608 :
1609 0 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1610 0 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix))
1611 0 : DO lmatrix = 1, qmatrix
1612 0 : DO kmatrix = 1, pmatrix
1613 0 : DO jmatrix = 1, mmatrix
1614 0 : DO imatrix = 1, nmatrix
1615 0 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1616 : END DO
1617 : END DO
1618 : END DO
1619 : END DO
1620 0 : END SUBROUTINE allocate_dbcsr_matrix_set_4d
1621 :
1622 : ! **************************************************************************************************
1623 : !> \brief Allocate and initialize a real matrix 5-dimensional set.
1624 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1625 : !> \param[in] nmatrix Size of set
1626 : !> \param mmatrix ...
1627 : !> \param pmatrix ...
1628 : !> \param qmatrix ...
1629 : !> \param smatrix ...
1630 : !> \par History
1631 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1632 : ! **************************************************************************************************
1633 0 : SUBROUTINE allocate_dbcsr_matrix_set_5d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix, smatrix)
1634 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), &
1635 : POINTER :: matrix_set
1636 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix, &
1637 : smatrix
1638 :
1639 : INTEGER :: hmatrix, imatrix, jmatrix, kmatrix, &
1640 : lmatrix
1641 :
1642 0 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1643 0 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix, smatrix))
1644 0 : DO hmatrix = 1, smatrix
1645 0 : DO lmatrix = 1, qmatrix
1646 0 : DO kmatrix = 1, pmatrix
1647 0 : DO jmatrix = 1, mmatrix
1648 0 : DO imatrix = 1, nmatrix
1649 0 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1650 : END DO
1651 : END DO
1652 : END DO
1653 : END DO
1654 : END DO
1655 0 : END SUBROUTINE allocate_dbcsr_matrix_set_5d
1656 :
1657 : ! **************************************************************************************************
1658 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1659 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1660 : !> \par History
1661 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1662 : ! **************************************************************************************************
1663 203308 : SUBROUTINE deallocate_dbcsr_matrix_set_1d(matrix_set)
1664 :
1665 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1666 :
1667 : INTEGER :: imatrix
1668 :
1669 203308 : IF (ASSOCIATED(matrix_set)) THEN
1670 1004464 : DO imatrix = 1, SIZE(matrix_set)
1671 1004464 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix)%matrix)
1672 : END DO
1673 201616 : DEALLOCATE (matrix_set)
1674 : END IF
1675 :
1676 203308 : END SUBROUTINE deallocate_dbcsr_matrix_set_1d
1677 :
1678 : ! **************************************************************************************************
1679 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1680 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1681 : !> \par History
1682 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1683 : ! **************************************************************************************************
1684 200225 : SUBROUTINE deallocate_dbcsr_matrix_set_2d(matrix_set)
1685 :
1686 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1687 :
1688 : INTEGER :: imatrix, jmatrix
1689 :
1690 200225 : IF (ASSOCIATED(matrix_set)) THEN
1691 2197320 : DO jmatrix = 1, SIZE(matrix_set, 2)
1692 4462695 : DO imatrix = 1, SIZE(matrix_set, 1)
1693 4265494 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix)%matrix)
1694 : END DO
1695 : END DO
1696 197201 : DEALLOCATE (matrix_set)
1697 : END IF
1698 200225 : END SUBROUTINE deallocate_dbcsr_matrix_set_2d
1699 :
1700 : ! **************************************************************************************************
1701 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1702 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1703 : !> \par History
1704 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1705 : ! **************************************************************************************************
1706 328 : SUBROUTINE deallocate_dbcsr_matrix_set_3d(matrix_set)
1707 :
1708 : TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1709 :
1710 : INTEGER :: imatrix, jmatrix, kmatrix
1711 :
1712 328 : IF (ASSOCIATED(matrix_set)) THEN
1713 668 : DO kmatrix = 1, SIZE(matrix_set, 3)
1714 2562 : DO jmatrix = 1, SIZE(matrix_set, 2)
1715 6022 : DO imatrix = 1, SIZE(matrix_set, 1)
1716 5682 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1717 : END DO
1718 : END DO
1719 : END DO
1720 328 : DEALLOCATE (matrix_set)
1721 : END IF
1722 328 : END SUBROUTINE deallocate_dbcsr_matrix_set_3d
1723 :
1724 : ! **************************************************************************************************
1725 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1726 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1727 : !> \par History
1728 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1729 : ! **************************************************************************************************
1730 0 : SUBROUTINE deallocate_dbcsr_matrix_set_4d(matrix_set)
1731 :
1732 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1733 :
1734 : INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1735 :
1736 0 : IF (ASSOCIATED(matrix_set)) THEN
1737 0 : DO lmatrix = 1, SIZE(matrix_set, 4)
1738 0 : DO kmatrix = 1, SIZE(matrix_set, 3)
1739 0 : DO jmatrix = 1, SIZE(matrix_set, 2)
1740 0 : DO imatrix = 1, SIZE(matrix_set, 1)
1741 0 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1742 : END DO
1743 : END DO
1744 : END DO
1745 : END DO
1746 0 : DEALLOCATE (matrix_set)
1747 : END IF
1748 0 : END SUBROUTINE deallocate_dbcsr_matrix_set_4d
1749 :
1750 : ! **************************************************************************************************
1751 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1752 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1753 : !> \par History
1754 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1755 : ! **************************************************************************************************
1756 0 : SUBROUTINE deallocate_dbcsr_matrix_set_5d(matrix_set)
1757 :
1758 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), &
1759 : POINTER :: matrix_set
1760 :
1761 : INTEGER :: hmatrix, imatrix, jmatrix, kmatrix, &
1762 : lmatrix
1763 :
1764 0 : IF (ASSOCIATED(matrix_set)) THEN
1765 0 : DO hmatrix = 1, SIZE(matrix_set, 5)
1766 0 : DO lmatrix = 1, SIZE(matrix_set, 4)
1767 0 : DO kmatrix = 1, SIZE(matrix_set, 3)
1768 0 : DO jmatrix = 1, SIZE(matrix_set, 2)
1769 0 : DO imatrix = 1, SIZE(matrix_set, 1)
1770 0 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1771 : END DO
1772 : END DO
1773 : END DO
1774 : END DO
1775 : END DO
1776 0 : DEALLOCATE (matrix_set)
1777 : END IF
1778 0 : END SUBROUTINE deallocate_dbcsr_matrix_set_5d
1779 :
1780 0 : END MODULE cp_dbcsr_operations
|
