Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 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 mathlib, ONLY: lcm, gcd
20 : USE cp_blacs_env, ONLY: cp_blacs_env_type
21 : USE cp_cfm_types, ONLY: cp_cfm_type
22 : USE dbcsr_api, ONLY: dbcsr_distribution_get, &
23 : dbcsr_convert_sizes_to_offsets, dbcsr_add, &
24 : dbcsr_complete_redistribute, dbcsr_copy, dbcsr_create, &
25 : dbcsr_deallocate_matrix, &
26 : dbcsr_desymmetrize, dbcsr_distribution_new, &
27 : dbcsr_get_data_type, dbcsr_get_info, dbcsr_get_matrix_type, &
28 : dbcsr_iterator_type, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
29 : dbcsr_iterator_start, dbcsr_iterator_stop, &
30 : dbcsr_multiply, dbcsr_norm, dbcsr_p_type, dbcsr_release, &
31 : dbcsr_reserve_all_blocks, dbcsr_scale, dbcsr_type, &
32 : dbcsr_valid_index, dbcsr_verify_matrix, &
33 : dbcsr_distribution_type, dbcsr_distribution_release, &
34 : dbcsr_norm_frobenius, &
35 : dbcsr_type_antisymmetric, dbcsr_type_complex_8, dbcsr_type_no_symmetry, dbcsr_type_real_8, &
36 : dbcsr_type_symmetric
37 : USE cp_fm_basic_linalg, ONLY: cp_fm_gemm
38 : USE cp_fm_struct, ONLY: cp_fm_struct_create, &
39 : cp_fm_struct_release, &
40 : cp_fm_struct_type
41 : USE cp_fm_types, ONLY: cp_fm_create, &
42 : cp_fm_get_info, &
43 : cp_fm_release, &
44 : cp_fm_to_fm, &
45 : cp_fm_type
46 : USE message_passing, ONLY: mp_para_env_type
47 : USE distribution_2d_types, ONLY: distribution_2d_get, &
48 : distribution_2d_type
49 : USE kinds, ONLY: dp, default_string_length
50 : USE message_passing, ONLY: mp_comm_type
51 :
52 : !$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num, omp_get_num_threads
53 : #include "base/base_uses.f90"
54 :
55 : IMPLICIT NONE
56 :
57 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_operations'
58 : LOGICAL, PARAMETER :: debug_mod = .FALSE.
59 :
60 : INTEGER, SAVE, PUBLIC :: max_elements_per_block = 32
61 :
62 : PUBLIC :: dbcsr_multiply_local
63 :
64 : ! CP2K API emulation
65 : PUBLIC :: copy_fm_to_dbcsr, copy_dbcsr_to_fm, &
66 : copy_dbcsr_to_cfm, copy_cfm_to_dbcsr, &
67 : cp_dbcsr_sm_fm_multiply, cp_dbcsr_plus_fm_fm_t, &
68 : copy_dbcsr_to_fm_bc, copy_fm_to_dbcsr_bc, cp_fm_to_dbcsr_row_template, &
69 : cp_dbcsr_m_by_n_from_template, cp_dbcsr_m_by_n_from_row_template, &
70 : dbcsr_create_dist_r_unrot
71 :
72 : ! distribution_2d_type compatibility
73 : PUBLIC :: cp_dbcsr_dist2d_to_dist
74 :
75 : PUBLIC :: dbcsr_copy_columns_hack
76 :
77 : ! matrix set
78 : PUBLIC :: dbcsr_allocate_matrix_set
79 : PUBLIC :: dbcsr_deallocate_matrix_set
80 :
81 : INTERFACE dbcsr_allocate_matrix_set
82 : #:for ii in range(1, 6)
83 : MODULE PROCEDURE allocate_dbcsr_matrix_set_${ii}$d
84 : #:endfor
85 : END INTERFACE
86 :
87 : INTERFACE dbcsr_deallocate_matrix_set
88 : #:for ii in range(1, 6)
89 : MODULE PROCEDURE deallocate_dbcsr_matrix_set_${ii}$d
90 : #:endfor
91 : END INTERFACE
92 :
93 : PRIVATE
94 :
95 : CONTAINS
96 :
97 : #:for fm, type, constr in [("fm", "REAL", "REAL"), ("cfm", "COMPLEX", "CMPLX")]
98 :
99 : ! **************************************************************************************************
100 : !> \brief Copy a BLACS matrix to a dbcsr matrix.
101 : !>
102 : !> real_matrix=beta*real_matrix+alpha*fm
103 : !> beta defaults to 0, alpha to 1
104 : !> \param[in] fm full matrix
105 : !> \param[out] matrix DBCSR matrix
106 : !> \param[in] keep_sparsity (optional) retains the sparsity of the input
107 : !> matrix
108 : !> \date 2009-10-13
109 : !> \par History
110 : !> 2009-10-13 rewritten based on copy_dbcsr_to_fm
111 : !> \author Urban Borstnik
112 : !> \version 2.0
113 : ! **************************************************************************************************
114 1048043 : SUBROUTINE copy_${fm}$_to_dbcsr(fm, matrix, keep_sparsity)
115 : TYPE(cp_${fm}$_type), INTENT(IN) :: fm
116 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix
117 : LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
118 :
119 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_${fm}$_to_dbcsr'
120 :
121 : TYPE(dbcsr_type) :: bc_mat
122 : INTEGER :: handle
123 :
124 1048043 : CALL timeset(routineN, handle)
125 :
126 1048043 : CALL copy_${fm}$_to_dbcsr_bc(fm, bc_mat)
127 1048043 : CALL dbcsr_complete_redistribute(bc_mat, matrix, keep_sparsity=keep_sparsity)
128 1048043 : CALL dbcsr_release(bc_mat)
129 :
130 1048043 : CALL timestop(handle)
131 1048043 : END SUBROUTINE copy_${fm}$_to_dbcsr
132 :
133 : ! **************************************************************************************************
134 : !> \brief Copy a BLACS matrix to a dbcsr matrix with a special block-cyclic distribution,
135 : !> which requires no complete redistribution.
136 : !> \param fm ...
137 : !> \param bc_mat ...
138 : ! **************************************************************************************************
139 1048413 : SUBROUTINE copy_${fm}$_to_dbcsr_bc(fm, bc_mat)
140 : TYPE(cp_${fm}$_type), INTENT(IN) :: fm
141 : TYPE(dbcsr_type), INTENT(INOUT) :: bc_mat
142 :
143 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_${fm}$_to_dbcsr_bc'
144 :
145 : INTEGER :: col, handle, ncol_block, ncol_global, nrow_block, nrow_global, row
146 1048413 : INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
147 1048413 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
148 1048413 : ${type}$ (KIND=dp), DIMENSION(:, :), POINTER :: fm_block, dbcsr_block
149 : TYPE(dbcsr_distribution_type) :: bc_dist
150 : TYPE(dbcsr_iterator_type) :: iter
151 1048413 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
152 :
153 1048413 : CALL timeset(routineN, handle)
154 :
155 : #:if (type=="REAL")
156 1045875 : IF (fm%use_sp) CPABORT("copy_${fm}$_to_dbcsr_bc: single precision not supported")
157 : #:endif
158 :
159 : ! Create processor grid
160 1048413 : pgrid => fm%matrix_struct%context%blacs2mpi
161 :
162 : ! Create a block-cyclic distribution compatible with the FM matrix.
163 1048413 : nrow_block = fm%matrix_struct%nrow_block
164 1048413 : ncol_block = fm%matrix_struct%ncol_block
165 1048413 : nrow_global = fm%matrix_struct%nrow_global
166 1048413 : ncol_global = fm%matrix_struct%ncol_global
167 1048413 : NULLIFY (col_blk_size, row_blk_size)
168 : CALL dbcsr_create_dist_block_cyclic(bc_dist, &
169 : nrows=nrow_global, ncolumns=ncol_global, & ! Actual full matrix size
170 : nrow_block=nrow_block, ncol_block=ncol_block, & ! BLACS parameters
171 : group_handle=fm%matrix_struct%para_env%get_handle(), pgrid=pgrid, &
172 1048413 : row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size) ! block-cyclic row/col sizes
173 :
174 : ! Create the block-cyclic DBCSR matrix
175 : CALL dbcsr_create(bc_mat, "Block-cyclic ", bc_dist, &
176 : dbcsr_type_no_symmetry, row_blk_size, col_blk_size, nze=0, &
177 1048413 : reuse_arrays=.TRUE., data_type=dbcsr_type_${type.lower()}$_8)
178 1048413 : CALL dbcsr_distribution_release(bc_dist)
179 :
180 : ! allocate all blocks
181 1048413 : CALL dbcsr_reserve_all_blocks(bc_mat)
182 :
183 1048413 : CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
184 :
185 : ! Copy the FM data to the block-cyclic DBCSR matrix. This step
186 : ! could be skipped with appropriate DBCSR index manipulation.
187 1048413 : fm_block => fm%local_data
188 : !$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
189 1048413 : !$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
190 : CALL dbcsr_iterator_start(iter, bc_mat)
191 : DO WHILE (dbcsr_iterator_blocks_left(iter))
192 : CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
193 : dbcsr_block(:, :) = fm_block(first_row(row):last_row(row), first_col(col):last_col(col))
194 : END DO
195 : CALL dbcsr_iterator_stop(iter)
196 : !$OMP END PARALLEL
197 :
198 1048413 : CALL timestop(handle)
199 2096826 : END SUBROUTINE copy_${fm}$_to_dbcsr_bc
200 :
201 : ! **************************************************************************************************
202 : !> \brief Copy a DBCSR matrix to a BLACS matrix
203 : !> \param[in] matrix DBCSR matrix
204 : !> \param[out] fm full matrix
205 : ! **************************************************************************************************
206 3899428 : SUBROUTINE copy_dbcsr_to_${fm}$ (matrix, fm)
207 : TYPE(dbcsr_type), INTENT(IN) :: matrix
208 : TYPE(cp_${fm}$_type), INTENT(IN) :: fm
209 :
210 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_dbcsr_to_${fm}$'
211 :
212 974857 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
213 : INTEGER :: handle, ncol_block, nfullcols_total, group_handle, &
214 : nfullrows_total, nrow_block
215 : TYPE(dbcsr_type) :: bc_mat, matrix_nosym
216 : TYPE(dbcsr_distribution_type) :: dist, bc_dist
217 : CHARACTER(len=default_string_length) :: name
218 974857 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
219 :
220 974857 : CALL timeset(routineN, handle)
221 :
222 : ! check compatibility
223 : CALL dbcsr_get_info(matrix, &
224 : name=name, &
225 : distribution=dist, &
226 : nfullrows_total=nfullrows_total, &
227 974857 : nfullcols_total=nfullcols_total)
228 :
229 974857 : CPASSERT(fm%matrix_struct%nrow_global == nfullrows_total)
230 974857 : CPASSERT(fm%matrix_struct%ncol_global == nfullcols_total)
231 :
232 : ! info about the full matrix
233 974857 : nrow_block = fm%matrix_struct%nrow_block
234 974857 : ncol_block = fm%matrix_struct%ncol_block
235 :
236 : ! Convert DBCSR to a block-cyclic
237 974857 : NULLIFY (col_blk_size, row_blk_size)
238 974857 : CALL dbcsr_distribution_get(dist, group=group_handle, pgrid=pgrid)
239 : CALL dbcsr_create_dist_block_cyclic(bc_dist, &
240 : nrows=nfullrows_total, ncolumns=nfullcols_total, &
241 : nrow_block=nrow_block, ncol_block=ncol_block, &
242 : group_handle=group_handle, pgrid=pgrid, &
243 974857 : row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
244 :
245 : CALL dbcsr_create(bc_mat, "Block-cyclic"//name, bc_dist, &
246 : dbcsr_type_no_symmetry, row_blk_size, col_blk_size, &
247 : nze=0, data_type=dbcsr_get_data_type(matrix), &
248 974857 : reuse_arrays=.TRUE.)
249 974857 : CALL dbcsr_distribution_release(bc_dist)
250 :
251 974857 : CALL dbcsr_create(matrix_nosym, template=matrix, matrix_type="N")
252 974857 : CALL dbcsr_desymmetrize(matrix, matrix_nosym)
253 974857 : CALL dbcsr_complete_redistribute(matrix_nosym, bc_mat)
254 974857 : CALL dbcsr_release(matrix_nosym)
255 :
256 974857 : CALL copy_dbcsr_to_${fm}$_bc(bc_mat, fm)
257 :
258 974857 : CALL dbcsr_release(bc_mat)
259 :
260 974857 : CALL timestop(handle)
261 974857 : END SUBROUTINE copy_dbcsr_to_${fm}$
262 :
263 : ! **************************************************************************************************
264 : !> \brief Copy a DBCSR_BLACS matrix to a BLACS matrix
265 : !> \param bc_mat DBCSR matrix
266 : !> \param[out] fm full matrix
267 : ! **************************************************************************************************
268 974857 : SUBROUTINE copy_dbcsr_to_${fm}$_bc(bc_mat, fm)
269 : TYPE(dbcsr_type), INTENT(IN) :: bc_mat
270 : TYPE(cp_${fm}$_type), INTENT(IN) :: fm
271 :
272 : CHARACTER(LEN=*), PARAMETER :: routineN = 'copy_dbcsr_to_${fm}$_bc'
273 :
274 : INTEGER :: col, handle, row
275 974857 : INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
276 974857 : ${type}$ (KIND=dp), DIMENSION(:, :), POINTER :: dbcsr_block, fm_block
277 : TYPE(dbcsr_iterator_type) :: iter
278 :
279 974857 : CALL timeset(routineN, handle)
280 :
281 : #:if (type=="REAL")
282 972319 : IF (fm%use_sp) CPABORT("copy_dbcsr_to_${fm}$_bc: single precision not supported")
283 : #:endif
284 :
285 974857 : CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
286 :
287 : ! Now copy data to the FM matrix
288 974857 : fm_block => fm%local_data
289 377007429 : fm_block = ${constr}$ (0.0, KIND=dp)
290 : !$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
291 974857 : !$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
292 : CALL dbcsr_iterator_start(iter, bc_mat)
293 : DO WHILE (dbcsr_iterator_blocks_left(iter))
294 : CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
295 : fm_block(first_row(row):last_row(row), first_col(col):last_col(col)) = dbcsr_block(:, :)
296 : END DO
297 : CALL dbcsr_iterator_stop(iter)
298 : !$OMP END PARALLEL
299 :
300 974857 : CALL timestop(handle)
301 1949714 : END SUBROUTINE copy_dbcsr_to_${fm}$_bc
302 :
303 : #:endfor
304 :
305 : ! **************************************************************************************************
306 : !> \brief Helper routine used to copy blocks from DBCSR into FM matrices and vice versa
307 : !> \param bc_mat ...
308 : !> \param first_row ...
309 : !> \param last_row ...
310 : !> \param first_col ...
311 : !> \param last_col ...
312 : !> \author Ole Schuett
313 : ! **************************************************************************************************
314 2023270 : SUBROUTINE calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
315 : TYPE(dbcsr_type), INTENT(IN) :: bc_mat
316 : INTEGER :: col, nblkcols_local, nblkcols_total, nblkrows_local, nblkrows_total, row
317 : INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, &
318 : last_row, local_col_sizes, &
319 : local_row_sizes
320 2023270 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, local_cols, local_rows, &
321 2023270 : row_blk_size
322 :
323 : CALL dbcsr_get_info(bc_mat, &
324 : nblkrows_total=nblkrows_total, &
325 : nblkcols_total=nblkcols_total, &
326 : nblkrows_local=nblkrows_local, &
327 : nblkcols_local=nblkcols_local, &
328 : local_rows=local_rows, &
329 : local_cols=local_cols, &
330 : row_blk_size=row_blk_size, &
331 2023270 : col_blk_size=col_blk_size)
332 :
333 : ! calculate first_row and last_row
334 6069298 : ALLOCATE (local_row_sizes(nblkrows_total))
335 6042632 : local_row_sizes(:) = 0
336 2023270 : IF (nblkrows_local .GE. 1) THEN
337 4051574 : DO row = 1, nblkrows_local
338 4051574 : local_row_sizes(local_rows(row)) = row_blk_size(local_rows(row))
339 : END DO
340 : END IF
341 10114302 : ALLOCATE (first_row(nblkrows_total), last_row(nblkrows_total))
342 2023270 : CALL dbcsr_convert_sizes_to_offsets(local_row_sizes, first_row, last_row)
343 :
344 : ! calculate first_col and last_col
345 6068028 : ALLOCATE (local_col_sizes(nblkcols_total))
346 4693444 : local_col_sizes(:) = 0
347 2023270 : IF (nblkcols_local .GE. 1) THEN
348 4691662 : DO col = 1, nblkcols_local
349 4691662 : local_col_sizes(local_cols(col)) = col_blk_size(local_cols(col))
350 : END DO
351 : END IF
352 10109222 : ALLOCATE (first_col(nblkcols_total), last_col(nblkcols_total))
353 2023270 : CALL dbcsr_convert_sizes_to_offsets(local_col_sizes, first_col, last_col)
354 :
355 4046540 : END SUBROUTINE calculate_fm_block_ranges
356 :
357 : ! **************************************************************************************************
358 : !> \brief hack for dbcsr_copy_columns
359 : !> \param matrix_b ...
360 : !> \param matrix_a ...
361 : !> \param ncol ...
362 : !> \param source_start ...
363 : !> \param target_start ...
364 : !> \param para_env ...
365 : !> \param blacs_env ...
366 : !> \author vw
367 : ! **************************************************************************************************
368 9048 : SUBROUTINE dbcsr_copy_columns_hack(matrix_b, matrix_a, &
369 : ncol, source_start, target_start, para_env, blacs_env)
370 :
371 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix_b
372 : TYPE(dbcsr_type), INTENT(IN) :: matrix_a
373 : INTEGER, INTENT(IN) :: ncol, source_start, target_start
374 : TYPE(mp_para_env_type), POINTER :: para_env
375 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
376 :
377 : INTEGER :: nfullcols_total, nfullrows_total
378 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
379 : TYPE(cp_fm_type) :: fm_matrix_a, fm_matrix_b
380 :
381 2262 : NULLIFY (fm_struct)
382 2262 : CALL dbcsr_get_info(matrix_a, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
383 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
384 2262 : ncol_global=nfullcols_total, para_env=para_env)
385 2262 : CALL cp_fm_create(fm_matrix_a, fm_struct, name="fm_matrix_a")
386 2262 : CALL cp_fm_struct_release(fm_struct)
387 :
388 2262 : CALL dbcsr_get_info(matrix_b, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
389 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
390 2262 : ncol_global=nfullcols_total, para_env=para_env)
391 2262 : CALL cp_fm_create(fm_matrix_b, fm_struct, name="fm_matrix_b")
392 2262 : CALL cp_fm_struct_release(fm_struct)
393 :
394 2262 : CALL copy_dbcsr_to_fm(matrix_a, fm_matrix_a)
395 2262 : CALL copy_dbcsr_to_fm(matrix_b, fm_matrix_b)
396 :
397 2262 : CALL cp_fm_to_fm(fm_matrix_a, fm_matrix_b, ncol, source_start, target_start)
398 :
399 2262 : CALL copy_fm_to_dbcsr(fm_matrix_b, matrix_b)
400 :
401 2262 : CALL cp_fm_release(fm_matrix_a)
402 2262 : CALL cp_fm_release(fm_matrix_b)
403 :
404 2262 : END SUBROUTINE dbcsr_copy_columns_hack
405 :
406 : ! **************************************************************************************************
407 : !> \brief Creates a DBCSR distribution from a distribution_2d
408 : !> \param[in] dist2d distribution_2d
409 : !> \param[out] dist DBCSR distribution
410 : !> \par History
411 : !> move form dbcsr_operation 01.2010
412 : ! **************************************************************************************************
413 8612 : SUBROUTINE cp_dbcsr_dist2d_to_dist(dist2d, dist)
414 : TYPE(distribution_2d_type), INTENT(IN), TARGET :: dist2d
415 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
416 :
417 8612 : INTEGER, DIMENSION(:, :), POINTER :: pgrid, col_dist_data, row_dist_data
418 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
419 : TYPE(mp_para_env_type), POINTER :: para_env
420 : TYPE(distribution_2d_type), POINTER :: dist2d_p
421 8612 : INTEGER, DIMENSION(:), POINTER :: row_dist, col_dist
422 :
423 8612 : dist2d_p => dist2d
424 : CALL distribution_2d_get(dist2d_p, &
425 : row_distribution=row_dist_data, &
426 : col_distribution=col_dist_data, &
427 8612 : blacs_env=blacs_env)
428 8612 : CALL blacs_env%get(para_env=para_env, blacs2mpi=pgrid)
429 :
430 : ! map to 1D arrays
431 8612 : row_dist => row_dist_data(:, 1)
432 8612 : col_dist => col_dist_data(:, 1)
433 : !row_cluster => row_dist_data(:, 2)
434 : !col_cluster => col_dist_data(:, 2)
435 :
436 : CALL dbcsr_distribution_new(dist, &
437 : group=para_env%get_handle(), pgrid=pgrid, &
438 : row_dist=row_dist, &
439 8612 : col_dist=col_dist)
440 :
441 8612 : END SUBROUTINE cp_dbcsr_dist2d_to_dist
442 :
443 : ! **************************************************************************************************
444 : !> \brief multiply a dbcsr with a replicated array
445 : !> c = alpha_scalar * A (dbscr) * b + c
446 : !> \param[in] matrix_a DBSCR matrxx
447 : !> \param[in] vec_b vectors b
448 : !> \param[inout] vec_c vectors c
449 : !> \param[in] ncol nbr of columns
450 : !> \param[in] alpha alpha
451 : !>
452 : ! **************************************************************************************************
453 0 : SUBROUTINE dbcsr_multiply_local(matrix_a, vec_b, vec_c, ncol, alpha)
454 : TYPE(dbcsr_type), INTENT(IN) :: matrix_a
455 : REAL(dp), DIMENSION(:, :), INTENT(IN) :: vec_b
456 : REAL(dp), DIMENSION(:, :), INTENT(INOUT) :: vec_c
457 : INTEGER, INTENT(in), OPTIONAL :: ncol
458 : REAL(dp), INTENT(IN), OPTIONAL :: alpha
459 :
460 : CHARACTER(LEN=*), PARAMETER :: routineN = 'dbcsr_multiply_local'
461 :
462 : INTEGER :: blk, col, coloff, my_ncol, row, rowoff, &
463 : timing_handle
464 : LOGICAL :: has_symm
465 : REAL(dp) :: my_alpha, my_alpha2
466 0 : REAL(dp), DIMENSION(:, :), POINTER :: data_d
467 : TYPE(dbcsr_iterator_type) :: iter
468 :
469 0 : CALL timeset(routineN, timing_handle)
470 :
471 0 : my_alpha = 1.0_dp
472 0 : IF (PRESENT(alpha)) my_alpha = alpha
473 :
474 0 : my_ncol = SIZE(vec_b, 2)
475 0 : IF (PRESENT(ncol)) my_ncol = ncol
476 :
477 0 : my_alpha2 = 0.0_dp
478 0 : IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric) my_alpha2 = my_alpha
479 0 : IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric) my_alpha2 = -my_alpha
480 :
481 : has_symm = (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric .OR. &
482 0 : dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric)
483 :
484 : !$OMP PARALLEL DEFAULT(NONE) SHARED(matrix_a,vec_b,vec_c,ncol,my_alpha2,my_alpha,my_ncol,has_symm) &
485 0 : !$OMP PRIVATE(iter,row,col,data_d,blk,rowoff,coloff)
486 : CALL dbcsr_iterator_start(iter, matrix_a, read_only=.TRUE., dynamic=.TRUE., dynamic_byrows=.TRUE.)
487 : DO WHILE (dbcsr_iterator_blocks_left(iter))
488 : CALL dbcsr_iterator_next_block(iter, row, col, data_d, blk, row_offset=rowoff, col_offset=coloff)
489 : IF (my_ncol .NE. 1) THEN
490 : CALL dgemm('N', 'N', &
491 : SIZE(data_d, 1), my_ncol, SIZE(data_d, 2), &
492 : my_alpha, data_d(1, 1), SIZE(data_d, 1), &
493 : vec_b(coloff, 1), SIZE(vec_b, 1), &
494 : 1.0_dp, vec_c(rowoff, 1), SIZE(vec_c, 1))
495 : ELSE
496 : CALL dgemv('N', SIZE(data_d, 1), SIZE(data_d, 2), &
497 : my_alpha, data_d(1, 1), SIZE(data_d, 1), &
498 : vec_b(coloff, 1), 1, &
499 : 1.0_dp, vec_c(rowoff, 1), 1)
500 : END IF
501 : END DO
502 : CALL dbcsr_iterator_stop(iter)
503 : !$OMP END PARALLEL
504 :
505 : ! FIXME ... in the symmetric case, the writes to vec_c depend on the column, not the row. This makes OMP-ing more difficult
506 : ! needs e.g. a buffer for vec_c and a reduction of that buffer.
507 0 : IF (has_symm) THEN
508 0 : CALL dbcsr_iterator_start(iter, matrix_a)
509 0 : DO WHILE (dbcsr_iterator_blocks_left(iter))
510 0 : CALL dbcsr_iterator_next_block(iter, row, col, data_d, blk, row_offset=rowoff, col_offset=coloff)
511 0 : IF (row .NE. col) THEN
512 0 : IF (my_ncol .NE. 1) THEN
513 : CALL dgemm('T', 'N', &
514 : SIZE(data_d, 2), my_ncol, SIZE(data_d, 1), &
515 : my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
516 : vec_b(rowoff, 1), SIZE(vec_b, 1), &
517 0 : 1.0_dp, vec_c(coloff, 1), SIZE(vec_c, 1))
518 : ELSE
519 : CALL dgemv('T', SIZE(data_d, 1), SIZE(data_d, 2), &
520 : my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
521 : vec_b(rowoff, 1), 1, &
522 0 : 1.0_dp, vec_c(coloff, 1), 1)
523 : END IF
524 : END IF
525 : END DO
526 0 : CALL dbcsr_iterator_stop(iter)
527 : END IF
528 :
529 0 : CALL timestop(timing_handle)
530 0 : END SUBROUTINE dbcsr_multiply_local
531 :
532 : ! **************************************************************************************************
533 : !> \brief multiply a dbcsr with a fm matrix
534 : !>
535 : !> For backwards compatibility with BLAS XGEMM, this routine supports
536 : !> the multiplication of matrices with incompatible dimensions.
537 : !>
538 : !> \param[in] matrix DBCSR matrix
539 : !> \param fm_in full matrix
540 : !> \param fm_out full matrix
541 : !> \param[in] ncol nbr of columns
542 : !> \param[in] alpha alpha
543 : !> \param[in] beta beta
544 : !>
545 : ! **************************************************************************************************
546 2692278 : SUBROUTINE cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
547 : TYPE(dbcsr_type), INTENT(IN) :: matrix
548 : TYPE(cp_fm_type), INTENT(IN) :: fm_in, fm_out
549 : INTEGER, INTENT(IN) :: ncol
550 : REAL(dp), INTENT(IN), OPTIONAL :: alpha, beta
551 :
552 : CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_dbcsr_sm_fm_multiply'
553 :
554 : INTEGER :: k_in, k_out, timing_handle, &
555 : timing_handle_mult, &
556 : a_ncol, a_nrow, b_ncol, b_nrow, c_ncol, c_nrow
557 448713 : INTEGER, DIMENSION(:), POINTER :: col_blk_size_right_in, &
558 448713 : col_blk_size_right_out, row_blk_size, &
559 : !row_cluster, col_cluster,&
560 448713 : row_dist, col_dist, col_blk_size
561 : TYPE(dbcsr_type) :: in, out
562 : REAL(dp) :: my_alpha, my_beta
563 : TYPE(dbcsr_distribution_type) :: dist, dist_right_in, product_dist
564 :
565 448713 : CALL timeset(routineN, timing_handle)
566 :
567 448713 : my_alpha = 1.0_dp
568 448713 : my_beta = 0.0_dp
569 448713 : IF (PRESENT(alpha)) my_alpha = alpha
570 448713 : IF (PRESENT(beta)) my_beta = beta
571 :
572 : ! TODO
573 448713 : CALL cp_fm_get_info(fm_in, ncol_global=b_ncol, nrow_global=b_nrow)
574 448713 : CALL cp_fm_get_info(fm_out, ncol_global=c_ncol, nrow_global=c_nrow)
575 448713 : CALL dbcsr_get_info(matrix, nfullrows_total=a_nrow, nfullcols_total=a_ncol)
576 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: A ", a_nrow, "x", a_ncol
577 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: B ", b_nrow, "x", b_ncol
578 : !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: C ", c_nrow, "x", c_ncol
579 :
580 448713 : CALL cp_fm_get_info(fm_out, ncol_global=k_out)
581 :
582 448713 : CALL cp_fm_get_info(fm_in, ncol_global=k_in)
583 : !write(*,*)routineN//" -----------------------------------"
584 : !IF (k_in .NE. k_out) &
585 : ! WRITE(*,'(3(A,I5,1X),2(A,F5.2,1X))')&
586 : ! routineN//" ncol", ncol,'k_in',k_in,'k_out',k_out,&
587 : ! 'alpha',my_alpha,'beta',my_beta
588 :
589 448713 : IF (ncol .GT. 0 .AND. k_out .GT. 0 .AND. k_in .GT. 0) THEN
590 447687 : CALL dbcsr_get_info(matrix, row_blk_size=row_blk_size, col_blk_size=col_blk_size, distribution=dist)
591 447687 : CALL dbcsr_create_dist_r_unrot(dist_right_in, dist, k_in, col_blk_size_right_in)
592 :
593 : CALL dbcsr_create(in, "D", dist_right_in, dbcsr_type_no_symmetry, &
594 447687 : col_blk_size, col_blk_size_right_in, nze=0)
595 :
596 447687 : CALL dbcsr_distribution_get(dist, row_dist=row_dist)
597 447687 : CALL dbcsr_distribution_get(dist_right_in, col_dist=col_dist)
598 : CALL dbcsr_distribution_new(product_dist, template=dist, &
599 447687 : row_dist=row_dist, col_dist=col_dist)
600 1343061 : ALLOCATE (col_blk_size_right_out(SIZE(col_blk_size_right_in)))
601 1817508 : col_blk_size_right_out = col_blk_size_right_in
602 447687 : CALL match_col_sizes(col_blk_size_right_out, col_blk_size_right_in, k_out)
603 :
604 : !if (k_in .ne. k_out) then
605 : ! write(*,*)routineN//" in cs", col_blk_size_right_in
606 : ! write(*,*)routineN//" out cs", col_blk_size_right_out
607 : !endif
608 :
609 : CALL dbcsr_create(out, "D", product_dist, dbcsr_type_no_symmetry, &
610 447687 : row_blk_size, col_blk_size_right_out, nze=0)
611 :
612 447687 : CALL copy_fm_to_dbcsr(fm_in, in)
613 447687 : IF (ncol .NE. k_out .OR. my_beta .NE. 0.0_dp) &
614 102998 : CALL copy_fm_to_dbcsr(fm_out, out)
615 :
616 447687 : CALL timeset(routineN//'_core', timing_handle_mult)
617 : CALL dbcsr_multiply("N", "N", my_alpha, matrix, in, my_beta, out, &
618 447687 : last_column=ncol)
619 447687 : CALL timestop(timing_handle_mult)
620 :
621 447687 : CALL copy_dbcsr_to_fm(out, fm_out)
622 :
623 447687 : CALL dbcsr_release(in)
624 447687 : CALL dbcsr_release(out)
625 447687 : DEALLOCATE (col_blk_size_right_in, col_blk_size_right_out)
626 447687 : CALL dbcsr_distribution_release(dist_right_in)
627 1790748 : CALL dbcsr_distribution_release(product_dist)
628 :
629 : END IF
630 :
631 448713 : CALL timestop(timing_handle)
632 :
633 448713 : END SUBROUTINE cp_dbcsr_sm_fm_multiply
634 :
635 : ! **************************************************************************************************
636 : !> \brief ...
637 : !> \param sizes1 ...
638 : !> \param sizes2 ...
639 : !> \param full_num ...
640 : ! **************************************************************************************************
641 447687 : SUBROUTINE match_col_sizes(sizes1, sizes2, full_num)
642 : INTEGER, DIMENSION(:), INTENT(INOUT) :: sizes1
643 : INTEGER, DIMENSION(:), INTENT(IN) :: sizes2
644 : INTEGER, INTENT(IN) :: full_num
645 :
646 : INTEGER :: left, n1, n2, p, rm, used
647 :
648 447687 : n1 = SIZE(sizes1)
649 447687 : n2 = SIZE(sizes2)
650 447687 : IF (n1 .NE. n2) &
651 0 : CPABORT("distributions must be equal!")
652 908754 : sizes1(1:n1) = sizes2(1:n1)
653 908754 : used = SUM(sizes1(1:n1))
654 : ! If sizes1 does not cover everything, then we increase the
655 : ! size of the last block; otherwise we reduce the blocks
656 : ! (from the end) until it is small enough.
657 447687 : IF (used .LT. full_num) THEN
658 0 : sizes1(n1) = sizes1(n1) + full_num - used
659 : ELSE
660 447687 : left = used - full_num
661 447687 : p = n1
662 447687 : DO WHILE (left .GT. 0 .AND. p .GT. 0)
663 0 : rm = MIN(left, sizes1(p))
664 0 : sizes1(p) = sizes1(p) - rm
665 0 : left = left - rm
666 0 : p = p - 1
667 : END DO
668 : END IF
669 447687 : END SUBROUTINE match_col_sizes
670 :
671 : ! **************************************************************************************************
672 : !> \brief performs the multiplication sparse_matrix+dense_mat*dens_mat^T
673 : !> if matrix_g is not explicitly given, matrix_v^T will be used
674 : !> this can be important to save the necessary redistribute for a
675 : !> different matrix_g and increase performance.
676 : !> \param sparse_matrix ...
677 : !> \param matrix_v ...
678 : !> \param matrix_g ...
679 : !> \param ncol ...
680 : !> \param alpha ...
681 : !> \param keep_sparsity Determines if the sparsity of sparse_matrix is retained
682 : !> by default it is TRUE
683 : !> \param symmetry_mode There are the following modes
684 : !> 1: sparse_matrix += 0.5*alpha*(v*g^T+g^T*v) (symmetric update)
685 : !> -1: sparse_matrix += 0.5*alpha*(v*g^T-g^T*v) (skewsymmetric update)
686 : !> else: sparse_matrix += alpha*v*g^T (no symmetry, default)
687 : !> saves some redistribution steps
688 : ! **************************************************************************************************
689 190016 : SUBROUTINE cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
690 : TYPE(dbcsr_type), INTENT(INOUT) :: sparse_matrix
691 : TYPE(cp_fm_type), INTENT(IN) :: matrix_v
692 : TYPE(cp_fm_type), OPTIONAL, INTENT(IN) :: matrix_g
693 : INTEGER, INTENT(IN) :: ncol
694 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: alpha
695 : LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
696 : INTEGER, INTENT(IN), OPTIONAL :: symmetry_mode
697 :
698 : CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_dbcsr_plus_fm_fm_t_native'
699 :
700 : INTEGER :: npcols, k, nao, timing_handle, data_type, my_symmetry_mode
701 190016 : INTEGER, DIMENSION(:), POINTER :: col_blk_size_left, &
702 190016 : col_dist_left, row_blk_size, row_dist
703 : LOGICAL :: check_product, my_keep_sparsity
704 : REAL(KIND=dp) :: my_alpha, norm
705 : TYPE(dbcsr_type) :: mat_g, mat_v, sparse_matrix2, &
706 : sparse_matrix3
707 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
708 : TYPE(cp_fm_type) :: fm_matrix
709 : TYPE(dbcsr_distribution_type) :: dist_left, sparse_dist
710 :
711 190016 : check_product = .FALSE.
712 :
713 190016 : CALL timeset(routineN, timing_handle)
714 :
715 190016 : my_keep_sparsity = .TRUE.
716 190016 : IF (PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
717 :
718 190016 : my_symmetry_mode = 0
719 190016 : IF (PRESENT(symmetry_mode)) my_symmetry_mode = symmetry_mode
720 :
721 190016 : NULLIFY (col_dist_left)
722 :
723 190016 : IF (ncol .GT. 0) THEN
724 188568 : IF (.NOT. dbcsr_valid_index(sparse_matrix)) &
725 0 : CPABORT("sparse_matrix must pre-exist")
726 : !
727 : ! Setup matrix_v
728 188568 : CALL cp_fm_get_info(matrix_v, ncol_global=k)
729 : !WRITE(*,*)routineN//'truncated mult k, ncol',k,ncol,' PRESENT (matrix_g)',PRESENT (matrix_g)
730 188568 : CALL dbcsr_get_info(sparse_matrix, distribution=sparse_dist)
731 188568 : CALL dbcsr_distribution_get(sparse_dist, npcols=npcols, row_dist=row_dist)
732 188568 : CALL create_bl_distribution(col_dist_left, col_blk_size_left, k, npcols)
733 : CALL dbcsr_distribution_new(dist_left, template=sparse_dist, &
734 188568 : row_dist=row_dist, col_dist=col_dist_left)
735 188568 : DEALLOCATE (col_dist_left)
736 188568 : CALL dbcsr_get_info(sparse_matrix, row_blk_size=row_blk_size, data_type=data_type)
737 : CALL dbcsr_create(mat_v, "DBCSR matrix_v", dist_left, dbcsr_type_no_symmetry, &
738 188568 : row_blk_size, col_blk_size_left, nze=0, data_type=data_type)
739 188568 : CALL copy_fm_to_dbcsr(matrix_v, mat_v)
740 188568 : CALL dbcsr_verify_matrix(mat_v)
741 : !
742 : ! Setup matrix_g
743 188568 : IF (PRESENT(matrix_g)) THEN
744 : CALL dbcsr_create(mat_g, "DBCSR matrix_g", dist_left, dbcsr_type_no_symmetry, &
745 86679 : row_blk_size, col_blk_size_left, data_type=data_type)
746 86679 : CALL copy_fm_to_dbcsr(matrix_g, mat_g)
747 : END IF
748 : !
749 188568 : DEALLOCATE (col_blk_size_left)
750 188568 : CALL dbcsr_distribution_release(dist_left)
751 : !
752 : !
753 : IF (check_product) THEN
754 : CALL cp_fm_get_info(matrix_v, nrow_global=nao)
755 : CALL cp_fm_struct_create(fm_struct_tmp, context=matrix_v%matrix_struct%context, nrow_global=nao, &
756 : ncol_global=nao, para_env=matrix_v%matrix_struct%para_env)
757 : CALL cp_fm_create(fm_matrix, fm_struct_tmp, name="fm matrix")
758 : CALL cp_fm_struct_release(fm_struct_tmp)
759 : CALL copy_dbcsr_to_fm(sparse_matrix, fm_matrix)
760 : CALL dbcsr_copy(sparse_matrix3, sparse_matrix)
761 : END IF
762 : !
763 188568 : my_alpha = 1.0_dp
764 188568 : IF (PRESENT(alpha)) my_alpha = alpha
765 188568 : IF (PRESENT(matrix_g)) THEN
766 86679 : IF (my_symmetry_mode == 1) THEN
767 : ! Symmetric mode
768 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
769 : 1.0_dp, sparse_matrix, &
770 : retain_sparsity=my_keep_sparsity, &
771 37232 : last_k=ncol)
772 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_g, mat_v, &
773 : 1.0_dp, sparse_matrix, &
774 : retain_sparsity=my_keep_sparsity, &
775 37232 : last_k=ncol)
776 49447 : ELSE IF (my_symmetry_mode == -1) THEN
777 : ! Skewsymmetric mode
778 : CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
779 : 1.0_dp, sparse_matrix, &
780 : retain_sparsity=my_keep_sparsity, &
781 2198 : last_k=ncol)
782 : CALL dbcsr_multiply("N", "T", -0.5_dp*my_alpha, mat_g, mat_v, &
783 : 1.0_dp, sparse_matrix, &
784 : retain_sparsity=my_keep_sparsity, &
785 2198 : last_k=ncol)
786 : ELSE
787 : ! Normal mode
788 : CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_g, &
789 : 1.0_dp, sparse_matrix, &
790 : retain_sparsity=my_keep_sparsity, &
791 47249 : last_k=ncol)
792 : END IF
793 : ELSE
794 : CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_v, &
795 : 1.0_dp, sparse_matrix, &
796 : retain_sparsity=my_keep_sparsity, &
797 101889 : last_k=ncol)
798 : END IF
799 :
800 : IF (check_product) THEN
801 : IF (PRESENT(matrix_g)) THEN
802 : IF (my_symmetry_mode == 1) THEN
803 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
804 : 1.0_dp, fm_matrix)
805 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_g, matrix_v, &
806 : 1.0_dp, fm_matrix)
807 : ELSE IF (my_symmetry_mode == -1) THEN
808 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
809 : 1.0_dp, fm_matrix)
810 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, -0.5_dp*my_alpha, matrix_g, matrix_v, &
811 : 1.0_dp, fm_matrix)
812 : ELSE
813 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_g, &
814 : 1.0_dp, fm_matrix)
815 : END IF
816 : ELSE
817 : CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_v, &
818 : 1.0_dp, fm_matrix)
819 : END IF
820 :
821 : CALL dbcsr_copy(sparse_matrix2, sparse_matrix)
822 : CALL dbcsr_scale(sparse_matrix2, alpha_scalar=0.0_dp)
823 : CALL copy_fm_to_dbcsr(fm_matrix, sparse_matrix2, keep_sparsity=my_keep_sparsity)
824 : CALL dbcsr_add(sparse_matrix2, sparse_matrix, alpha_scalar=1.0_dp, &
825 : beta_scalar=-1.0_dp)
826 : CALL dbcsr_norm(sparse_matrix2, which_norm=dbcsr_norm_frobenius, &
827 : norm_scalar=norm)
828 : WRITE (*, *) 'nao=', nao, ' k=', k, ' ncol=', ncol, ' my_alpha=', my_alpha
829 : WRITE (*, *) 'PRESENT (matrix_g)', PRESENT(matrix_g)
830 : WRITE (*, *) 'matrix_type=', dbcsr_get_matrix_type(sparse_matrix)
831 : WRITE (*, *) 'norm(sm+alpha*v*g^t - fm+alpha*v*g^t)/n=', norm/REAL(nao, dp)
832 : IF (norm/REAL(nao, dp) .GT. 1e-12_dp) THEN
833 : !WRITE(*,*) 'fm_matrix'
834 : !DO j=1,SIZE(fm_matrix%local_data,2)
835 : ! DO i=1,SIZE(fm_matrix%local_data,1)
836 : ! WRITE(*,'(A,I3,A,I3,A,E26.16,A)') 'a(',i,',',j,')=',fm_matrix%local_data(i,j),';'
837 : ! ENDDO
838 : !ENDDO
839 : !WRITE(*,*) 'mat_v'
840 : !CALL dbcsr_print(mat_v,matlab_format=.TRUE.)
841 : !WRITE(*,*) 'mat_g'
842 : !CALL dbcsr_print(mat_g,matlab_format=.TRUE.)
843 : !WRITE(*,*) 'sparse_matrix'
844 : !CALL dbcsr_print(sparse_matrix,matlab_format=.TRUE.)
845 : !WRITE(*,*) 'sparse_matrix2 (-sm + sparse(fm))'
846 : !CALL dbcsr_print(sparse_matrix2,matlab_format=.TRUE.)
847 : !WRITE(*,*) 'sparse_matrix3 (copy of sm input)'
848 : !CALL dbcsr_print(sparse_matrix3,matlab_format=.TRUE.)
849 : !stop
850 : END IF
851 : CALL dbcsr_release(sparse_matrix2)
852 : CALL dbcsr_release(sparse_matrix3)
853 : CALL cp_fm_release(fm_matrix)
854 : END IF
855 188568 : CALL dbcsr_release(mat_v)
856 188568 : IF (PRESENT(matrix_g)) CALL dbcsr_release(mat_g)
857 : END IF
858 190016 : CALL timestop(timing_handle)
859 :
860 190016 : END SUBROUTINE cp_dbcsr_plus_fm_fm_t
861 :
862 : ! **************************************************************************************************
863 : !> \brief Utility function to copy a specially shaped fm to dbcsr_matrix
864 : !> The result matrix will be the matrix in dbcsr format
865 : !> with the row blocks sizes according to the block_sizes of the template
866 : !> and the col blocks sizes evenly blocked with the internal dbcsr conversion
867 : !> size (32 is the current default)
868 : !> \param matrix ...
869 : !> \param fm_in ...
870 : !> \param template ...
871 : ! **************************************************************************************************
872 19784 : SUBROUTINE cp_fm_to_dbcsr_row_template(matrix, fm_in, template)
873 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix
874 : TYPE(cp_fm_type), INTENT(IN) :: fm_in
875 : TYPE(dbcsr_type), INTENT(IN) :: template
876 :
877 : INTEGER :: k_in, data_type
878 4946 : INTEGER, DIMENSION(:), POINTER :: col_blk_size_right_in, row_blk_size
879 : TYPE(dbcsr_distribution_type) :: tmpl_dist, dist_right_in
880 :
881 4946 : CALL cp_fm_get_info(fm_in, ncol_global=k_in)
882 :
883 4946 : CALL dbcsr_get_info(template, distribution=tmpl_dist)
884 4946 : CALL dbcsr_create_dist_r_unrot(dist_right_in, tmpl_dist, k_in, col_blk_size_right_in)
885 4946 : CALL dbcsr_get_info(template, row_blk_size=row_blk_size, data_type=data_type)
886 : CALL dbcsr_create(matrix, "D", dist_right_in, dbcsr_type_no_symmetry, &
887 4946 : row_blk_size, col_blk_size_right_in, nze=0, data_type=data_type)
888 :
889 4946 : CALL copy_fm_to_dbcsr(fm_in, matrix)
890 4946 : DEALLOCATE (col_blk_size_right_in)
891 4946 : CALL dbcsr_distribution_release(dist_right_in)
892 :
893 4946 : END SUBROUTINE cp_fm_to_dbcsr_row_template
894 :
895 : ! **************************************************************************************************
896 : !> \brief Utility function to create an arbitrary shaped dbcsr matrix
897 : !> with the same processor grid as the template matrix
898 : !> both row sizes and col sizes are evenly blocked with the internal
899 : !> dbcsr_conversion size (32 is the current default)
900 : !> \param matrix dbcsr matrix to be created
901 : !> \param template template dbcsr matrix giving its mp_env
902 : !> \param m global row size of output matrix
903 : !> \param n global col size of output matrix
904 : !> \param sym ...
905 : !> \param data_type ...
906 : ! **************************************************************************************************
907 285136 : SUBROUTINE cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym, data_type)
908 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
909 : INTEGER, INTENT(IN) :: m, n
910 : CHARACTER, OPTIONAL, INTENT(IN) :: sym
911 : INTEGER, OPTIONAL, INTENT(IN) :: data_type
912 :
913 : CHARACTER :: mysym
914 : INTEGER :: my_data_type, nprows, npcols
915 142568 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, &
916 142568 : col_dist, row_blk_size, &
917 142568 : row_dist
918 : TYPE(dbcsr_distribution_type) :: tmpl_dist, dist_m_n
919 :
920 : CALL dbcsr_get_info(template, &
921 : matrix_type=mysym, &
922 : data_type=my_data_type, &
923 142568 : distribution=tmpl_dist)
924 :
925 142568 : IF (PRESENT(sym)) mysym = sym
926 142568 : IF (PRESENT(data_type)) my_data_type = data_type
927 :
928 142568 : NULLIFY (row_dist, col_dist)
929 142568 : NULLIFY (row_blk_size, col_blk_size)
930 : !NULLIFY (row_cluster, col_cluster)
931 :
932 142568 : CALL dbcsr_distribution_get(tmpl_dist, nprows=nprows, npcols=npcols)
933 142568 : CALL create_bl_distribution(row_dist, row_blk_size, m, nprows)
934 142568 : CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
935 : CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
936 : row_dist=row_dist, col_dist=col_dist, &
937 : !row_cluster=row_cluster, col_cluster=col_cluster, &
938 142568 : reuse_arrays=.TRUE.)
939 :
940 : CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, &
941 : row_blk_size, col_blk_size, nze=0, data_type=my_data_type, &
942 142568 : reuse_arrays=.TRUE.)
943 142568 : CALL dbcsr_distribution_release(dist_m_n)
944 :
945 142568 : END SUBROUTINE cp_dbcsr_m_by_n_from_template
946 :
947 : ! **************************************************************************************************
948 : !> \brief Utility function to create dbcsr matrix, m x n matrix (n arbitrary)
949 : !> with the same processor grid and row distribution as the template matrix
950 : !> col sizes are evenly blocked with the internal
951 : !> dbcsr_conversion size (32 is the current default)
952 : !> \param matrix dbcsr matrix to be created
953 : !> \param template template dbcsr matrix giving its mp_env
954 : !> \param n global col size of output matrix
955 : !> \param sym ...
956 : !> \param data_type ...
957 : ! **************************************************************************************************
958 504948 : SUBROUTINE cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym, data_type)
959 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
960 : INTEGER :: n
961 : CHARACTER, OPTIONAL :: sym
962 : INTEGER, OPTIONAL :: data_type
963 :
964 : CHARACTER :: mysym
965 : INTEGER :: my_data_type, npcols
966 126237 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
967 126237 : row_dist
968 : TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
969 :
970 252474 : mysym = dbcsr_get_matrix_type(template)
971 126237 : IF (PRESENT(sym)) mysym = sym
972 126237 : my_data_type = dbcsr_get_data_type(template)
973 126237 : IF (PRESENT(data_type)) my_data_type = data_type
974 :
975 126237 : CALL dbcsr_get_info(template, distribution=tmpl_dist)
976 : CALL dbcsr_distribution_get(tmpl_dist, &
977 : npcols=npcols, &
978 126237 : row_dist=row_dist)
979 :
980 126237 : NULLIFY (col_dist, col_blk_size)
981 126237 : CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
982 : CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
983 126237 : row_dist=row_dist, col_dist=col_dist)
984 :
985 126237 : CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
986 : CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, &
987 126237 : row_blk_size, col_blk_size, nze=0, data_type=my_data_type)
988 :
989 126237 : DEALLOCATE (col_dist, col_blk_size)
990 126237 : CALL dbcsr_distribution_release(dist_m_n)
991 :
992 126237 : END SUBROUTINE cp_dbcsr_m_by_n_from_row_template
993 :
994 : ! **************************************************************************************************
995 : !> \brief Distributes elements into blocks and into bins
996 : !>
997 : !> \param[out] block_distribution block distribution to bins
998 : !> \param[out] block_size sizes of blocks
999 : !> \param[in] nelements number of elements to bin
1000 : !> \param[in] nbins number of bins
1001 : !> \par Term clarification
1002 : !> An example: blocks are atom blocks and bins are process rows/columns.
1003 : ! **************************************************************************************************
1004 1052574 : SUBROUTINE create_bl_distribution(block_distribution, &
1005 : block_size, nelements, nbins)
1006 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: block_distribution, block_size
1007 : INTEGER, INTENT(IN) :: nelements, nbins
1008 :
1009 : CHARACTER(len=*), PARAMETER :: routineN = 'create_bl_distribution', &
1010 : routineP = moduleN//':'//routineN
1011 :
1012 : INTEGER :: bin, blk_layer, element_stack, els, &
1013 : estimated_blocks, max_blocks_per_bin, &
1014 : nblks, nblocks, stat
1015 1052574 : INTEGER, DIMENSION(:), POINTER :: blk_dist, blk_sizes
1016 :
1017 : ! ---------------------------------------------------------------------------
1018 :
1019 1052574 : NULLIFY (block_distribution)
1020 1052574 : NULLIFY (block_size)
1021 : ! Define the sizes on which we build the distribution.
1022 1052574 : IF (nelements .GT. 0) THEN
1023 :
1024 1043864 : nblocks = CEILING(REAL(nelements, KIND=dp)/REAL(max_elements_per_block, KIND=dp))
1025 1043864 : max_blocks_per_bin = CEILING(REAL(nblocks, KIND=dp)/REAL(nbins, KIND=dp))
1026 :
1027 : IF (debug_mod) THEN
1028 : WRITE (*, '(1X,A,1X,A,I7,A,I7,A)') routineP, "For", nelements, &
1029 : " elements and", nbins, " bins"
1030 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1031 : max_elements_per_block, " max elements per block"
1032 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1033 : nblocks, " blocks"
1034 : WRITE (*, '(1X,A,1X,A,I7,A)') routineP, "There are", &
1035 : max_blocks_per_bin, " max blocks/bin"
1036 : END IF
1037 :
1038 1043864 : estimated_blocks = max_blocks_per_bin*nbins
1039 3131592 : ALLOCATE (blk_dist(estimated_blocks), stat=stat)
1040 1043864 : IF (stat /= 0) &
1041 0 : CPABORT("blk_dist")
1042 2087728 : ALLOCATE (blk_sizes(estimated_blocks), stat=stat)
1043 : IF (stat /= 0) &
1044 0 : CPABORT("blk_sizes")
1045 1043864 : element_stack = 0
1046 1043864 : nblks = 0
1047 2148116 : DO blk_layer = 1, max_blocks_per_bin
1048 3376304 : DO bin = 0, nbins - 1
1049 1228188 : els = MIN(max_elements_per_block, nelements - element_stack)
1050 2332440 : IF (els .GT. 0) THEN
1051 1115668 : element_stack = element_stack + els
1052 1115668 : nblks = nblks + 1
1053 1115668 : blk_dist(nblks) = bin
1054 1115668 : blk_sizes(nblks) = els
1055 : IF (debug_mod) WRITE (*, '(1X,A,I5,A,I5,A,I5)') routineP//" Assigning", &
1056 : els, " elements as block", nblks, " to bin", bin
1057 : END IF
1058 : END DO
1059 : END DO
1060 : ! Create the output arrays.
1061 1043864 : IF (nblks .EQ. estimated_blocks) THEN
1062 931344 : block_distribution => blk_dist
1063 931344 : block_size => blk_sizes
1064 : ELSE
1065 337560 : ALLOCATE (block_distribution(nblks), stat=stat)
1066 : IF (stat /= 0) &
1067 0 : CPABORT("blk_dist")
1068 455600 : block_distribution(:) = blk_dist(1:nblks)
1069 112520 : DEALLOCATE (blk_dist)
1070 225040 : ALLOCATE (block_size(nblks), stat=stat)
1071 : IF (stat /= 0) &
1072 0 : CPABORT("blk_sizes")
1073 455600 : block_size(:) = blk_sizes(1:nblks)
1074 112520 : DEALLOCATE (blk_sizes)
1075 : END IF
1076 : ELSE
1077 8710 : ALLOCATE (block_distribution(0), stat=stat)
1078 : IF (stat /= 0) &
1079 0 : CPABORT("blk_dist")
1080 8710 : ALLOCATE (block_size(0), stat=stat)
1081 : IF (stat /= 0) &
1082 0 : CPABORT("blk_sizes")
1083 : END IF
1084 : 1579 FORMAT(I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5, 1X, I5)
1085 : IF (debug_mod) THEN
1086 : WRITE (*, '(1X,A,A)') routineP//" Distribution"
1087 : WRITE (*, 1579) block_distribution(:)
1088 : WRITE (*, '(1X,A,A)') routineP//" Sizes"
1089 : WRITE (*, 1579) block_size(:)
1090 : END IF
1091 1052574 : END SUBROUTINE create_bl_distribution
1092 :
1093 : ! **************************************************************************************************
1094 : !> \brief Creates a new distribution for the right matrix in a matrix
1095 : !> multiplication with unrotated grid.
1096 : !> \param[out] dist_right new distribution for the right matrix
1097 : !> \param[in] dist_left the distribution of the left matrix
1098 : !> \param[in] ncolumns number of columns in right matrix
1099 : !> \param[out] right_col_blk_sizes sizes of blocks in the created column
1100 : !> \par The new row distribution for the right matrix is the same as the row
1101 : !> distribution of the left matrix, while the column distribution is
1102 : !> created so that it is appropriate to the parallel environment.
1103 : ! **************************************************************************************************
1104 452633 : SUBROUTINE dbcsr_create_dist_r_unrot(dist_right, dist_left, ncolumns, &
1105 : right_col_blk_sizes)
1106 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist_right
1107 : TYPE(dbcsr_distribution_type), INTENT(IN) :: dist_left
1108 : INTEGER, INTENT(IN) :: ncolumns
1109 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: right_col_blk_sizes
1110 :
1111 : INTEGER :: multiplicity, nimages, ncols, nprows, npcols
1112 : INTEGER, ALLOCATABLE, DIMENSION(:) :: tmp_images
1113 452633 : INTEGER, DIMENSION(:), POINTER :: right_col_dist, right_row_dist, old_col_dist, dummy
1114 :
1115 : CALL dbcsr_distribution_get(dist_left, &
1116 : ncols=ncols, &
1117 : col_dist=old_col_dist, &
1118 : nprows=nprows, &
1119 452633 : npcols=npcols)
1120 :
1121 : ! Create the column distribution
1122 452633 : CALL create_bl_distribution(right_col_dist, right_col_blk_sizes, ncolumns, npcols)
1123 : ! Create an even row distribution.
1124 1810532 : ALLOCATE (right_row_dist(ncols), tmp_images(ncols))
1125 452633 : nimages = lcm(nprows, npcols)/nprows
1126 452633 : multiplicity = nprows/gcd(nprows, npcols)
1127 452633 : CALL rebin_distribution(right_row_dist, tmp_images, old_col_dist, nprows, multiplicity, nimages)
1128 :
1129 452633 : NULLIFY (dummy)
1130 : CALL dbcsr_distribution_new(dist_right, &
1131 : template=dist_left, &
1132 : row_dist=right_row_dist, &
1133 : col_dist=right_col_dist, &
1134 : !row_cluster=dummy,&
1135 : !col_cluster=dummy,&
1136 452633 : reuse_arrays=.TRUE.)
1137 452633 : DEALLOCATE (tmp_images)
1138 452633 : END SUBROUTINE dbcsr_create_dist_r_unrot
1139 :
1140 : ! **************************************************************************************************
1141 : !> \brief Makes new distribution with decimation and multiplicity
1142 : !> \param[out] new_bins new real distribution
1143 : !> \param[out] images new image distribution
1144 : !> \param[in] source_bins Basis for the new distribution and images
1145 : !> \param[in] nbins number of bins in the new real distribution
1146 : !> \param[in] multiplicity multiplicity
1147 : !> \param[in] nimages number of images in the new distribution
1148 : !> \par Definition of multiplicity and nimages
1149 : !> Multiplicity and decimation (number of images) are used to
1150 : !> match process grid coordinates on non-square process
1151 : !> grids. Given source_nbins and target_nbins, their relation is
1152 : !> source_nbins * target_multiplicity
1153 : !> = target_nbins * target_nimages.
1154 : !> It is best when both multiplicity and nimages are small. To
1155 : !> get these two factors, then, one can use the following formulas:
1156 : !> nimages = lcm(source_nbins, target_nbins) / target_nbins
1157 : !> multiplicity = target_nbins / gcd(source_nbins, target_nbins)
1158 : !> from the target's point of view (nimages = target_nimages).
1159 : !> \par Mapping
1160 : !> The new distribution comprises of real bins and images within
1161 : !> bins. These can be view as target_nbins*nimages virtual
1162 : !> columns. These same virtual columns are also
1163 : !> source_nbins*multiplicity in number. Therefore these virtual
1164 : !> columns are mapped from source_nbins*multiplicity onto
1165 : !> target_bins*nimages (each target bin has nimages images):
1166 : !> Source 4: |1 2 3|4 5 6|7 8 9|A B C| (4*3)
1167 : !> Target 6: |1 2|3 4|5 6|7 8|9 A|B C| (6*2)
1168 : !> multiplicity=3, nimages=2, 12 virtual columns (1-C).
1169 : !> Source bin elements are evenly mapped into one of multiplicity
1170 : !> virtual columns. Other (non-even, block-size aware) mappings
1171 : !> could be better.
1172 : ! **************************************************************************************************
1173 452633 : SUBROUTINE rebin_distribution(new_bins, images, source_bins, &
1174 : nbins, multiplicity, nimages)
1175 : INTEGER, DIMENSION(:), INTENT(OUT) :: new_bins, images
1176 : INTEGER, DIMENSION(:), INTENT(IN) :: source_bins
1177 : INTEGER, INTENT(IN) :: nbins, multiplicity, nimages
1178 :
1179 : INTEGER :: bin, i, old_nbins, virtual_bin
1180 452633 : INTEGER, ALLOCATABLE, DIMENSION(:) :: bin_multiplier
1181 :
1182 : ! ---------------------------------------------------------------------------
1183 :
1184 452633 : IF (MOD(nbins*nimages, multiplicity) .NE. 0) &
1185 0 : CPWARN("mulitplicity is not divisor of new process grid coordinate")
1186 452633 : old_nbins = (nbins*nimages)/multiplicity
1187 1357899 : ALLOCATE (bin_multiplier(0:old_nbins - 1))
1188 905266 : bin_multiplier(:) = 0
1189 2148769 : DO i = 1, SIZE(new_bins)
1190 1696136 : IF (i .LE. SIZE(source_bins)) THEN
1191 1696136 : bin = source_bins(i)
1192 : ELSE
1193 : ! Fill remainder with a cyclic distribution
1194 0 : bin = MOD(i, old_nbins)
1195 : END IF
1196 1696136 : virtual_bin = bin*multiplicity + bin_multiplier(bin)
1197 1696136 : new_bins(i) = virtual_bin/nimages
1198 1696136 : images(i) = 1 + MOD(virtual_bin, nimages)
1199 1696136 : bin_multiplier(bin) = bin_multiplier(bin) + 1
1200 2148769 : IF (bin_multiplier(bin) .GE. multiplicity) THEN
1201 724136 : bin_multiplier(bin) = 0
1202 : END IF
1203 : END DO
1204 452633 : END SUBROUTINE rebin_distribution
1205 :
1206 : ! **************************************************************************************************
1207 : !> \brief Creates a block-cyclic compatible distribution
1208 : !>
1209 : !> All blocks in a dimension, except for possibly the last
1210 : !> block, have the same size.
1211 : !> \param[out] dist the elemental distribution
1212 : !> \param[in] nrows number of full rows
1213 : !> \param[in] ncolumns number of full columns
1214 : !> \param[in] nrow_block size of row blocks
1215 : !> \param[in] ncol_block size of column blocks
1216 : !> \param[in] mp_env multiprocess environment
1217 : !> \param[out] row_blk_sizes row block sizes
1218 : !> \param[out] col_blk_sizes column block sizes
1219 : ! **************************************************************************************************
1220 2023270 : SUBROUTINE dbcsr_create_dist_block_cyclic(dist, nrows, ncolumns, &
1221 : nrow_block, ncol_block, group_handle, pgrid, row_blk_sizes, col_blk_sizes)
1222 : TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
1223 : INTEGER, INTENT(IN) :: nrows, ncolumns, nrow_block, ncol_block
1224 : INTEGER, INTENT(IN) :: group_handle
1225 : INTEGER, DIMENSION(:, :), POINTER :: pgrid
1226 : INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: row_blk_sizes, col_blk_sizes
1227 :
1228 : CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_create_dist_block_cyclic'
1229 :
1230 : INTEGER :: nblkcols, nblkrows, npcols, nprows, &
1231 : pdim, sz
1232 2023270 : INTEGER, DIMENSION(:), POINTER :: cd_data, rd_data
1233 :
1234 : ! Row sizes
1235 2023270 : IF (nrow_block .EQ. 0) THEN
1236 : nblkrows = 0
1237 : sz = 0
1238 : ELSE
1239 2023270 : nblkrows = nrows/nrow_block
1240 2023270 : sz = MOD(nrows, nrow_block)
1241 : END IF
1242 2023270 : IF (sz .GT. 0) nblkrows = nblkrows + 1
1243 10114814 : ALLOCATE (row_blk_sizes(nblkrows), rd_data(nblkrows))
1244 6042632 : row_blk_sizes = nrow_block
1245 2023270 : IF (sz .NE. 0) row_blk_sizes(nblkrows) = sz
1246 :
1247 : ! Column sizes
1248 2023270 : IF (ncol_block .EQ. 0) THEN
1249 : nblkcols = 0
1250 : sz = 0
1251 : ELSE
1252 2023270 : nblkcols = ncolumns/ncol_block
1253 2023270 : sz = MOD(ncolumns, ncol_block)
1254 : END IF
1255 2023270 : IF (sz .GT. 0) nblkcols = nblkcols + 1
1256 10111004 : ALLOCATE (col_blk_sizes(nblkcols), cd_data(nblkcols))
1257 4693444 : col_blk_sizes = ncol_block
1258 2023270 : IF (sz .NE. 0) col_blk_sizes(nblkcols) = sz
1259 : !
1260 : IF (debug_mod) THEN
1261 : WRITE (*, *) routineN//" nrows,nrow_block,nblkrows=", &
1262 : nrows, nrow_block, nblkrows
1263 : WRITE (*, *) routineN//" ncols,ncol_block,nblkcols=", &
1264 : ncolumns, ncol_block, nblkcols
1265 : END IF
1266 : ! Calculate process row distribution
1267 2023270 : nprows = SIZE(pgrid, 1)
1268 6004324 : DO pdim = 0, MIN(nprows - 1, nblkrows - 1)
1269 10023686 : rd_data(1 + pdim:nblkrows:nprows) = pdim
1270 : END DO
1271 : ! Calculate process column distribution
1272 2023270 : npcols = SIZE(pgrid, 2)
1273 4044758 : DO pdim = 0, MIN(npcols - 1, nblkcols - 1)
1274 6714932 : cd_data(1 + pdim:nblkcols:npcols) = pdim
1275 : END DO
1276 : !
1277 : IF (debug_mod) THEN
1278 : WRITE (*, *) routineN//" row_dist", &
1279 : rd_data
1280 : WRITE (*, *) routineN//" col_dist", &
1281 : cd_data
1282 : END IF
1283 : !
1284 : CALL dbcsr_distribution_new(dist, &
1285 : group=group_handle, pgrid=pgrid, &
1286 : row_dist=rd_data, &
1287 : col_dist=cd_data, &
1288 2023270 : reuse_arrays=.TRUE.)
1289 :
1290 2023270 : END SUBROUTINE dbcsr_create_dist_block_cyclic
1291 :
1292 : ! **************************************************************************************************
1293 : !> \brief Allocate and initialize a real matrix 1-dimensional set.
1294 : !> \param[in,out] matrix_set Set containing the DBCSR matrices
1295 : !> \param[in] nmatrix Size of set
1296 : !> \par History
1297 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1298 : ! **************************************************************************************************
1299 179989 : SUBROUTINE allocate_dbcsr_matrix_set_1d(matrix_set, nmatrix)
1300 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1301 : INTEGER, INTENT(IN) :: nmatrix
1302 :
1303 : INTEGER :: imatrix
1304 :
1305 179989 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1306 843656 : ALLOCATE (matrix_set(nmatrix))
1307 483678 : DO imatrix = 1, nmatrix
1308 483678 : NULLIFY (matrix_set(imatrix)%matrix)
1309 : END DO
1310 179989 : END SUBROUTINE allocate_dbcsr_matrix_set_1d
1311 :
1312 : ! **************************************************************************************************
1313 : !> \brief Allocate and initialize a real matrix 2-dimensional set.
1314 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1315 : !> \param[in] nmatrix Size of set
1316 : !> \param mmatrix ...
1317 : !> \par History
1318 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1319 : ! **************************************************************************************************
1320 134858 : SUBROUTINE allocate_dbcsr_matrix_set_2d(matrix_set, nmatrix, mmatrix)
1321 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1322 : INTEGER, INTENT(IN) :: nmatrix, mmatrix
1323 :
1324 : INTEGER :: imatrix, jmatrix
1325 :
1326 134858 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1327 1846856 : ALLOCATE (matrix_set(nmatrix, mmatrix))
1328 735438 : DO jmatrix = 1, mmatrix
1329 1442282 : DO imatrix = 1, nmatrix
1330 1307424 : NULLIFY (matrix_set(imatrix, jmatrix)%matrix)
1331 : END DO
1332 : END DO
1333 134858 : END SUBROUTINE allocate_dbcsr_matrix_set_2d
1334 :
1335 : ! **************************************************************************************************
1336 : !> \brief Allocate and initialize a real matrix 3-dimensional set.
1337 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1338 : !> \param[in] nmatrix Size of set
1339 : !> \param mmatrix ...
1340 : !> \param pmatrix ...
1341 : !> \par History
1342 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1343 : ! **************************************************************************************************
1344 0 : SUBROUTINE allocate_dbcsr_matrix_set_3d(matrix_set, nmatrix, mmatrix, pmatrix)
1345 : TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1346 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix
1347 :
1348 : INTEGER :: imatrix, jmatrix, kmatrix
1349 :
1350 0 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1351 0 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix))
1352 0 : DO kmatrix = 1, pmatrix
1353 0 : DO jmatrix = 1, mmatrix
1354 0 : DO imatrix = 1, nmatrix
1355 0 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1356 : END DO
1357 : END DO
1358 : END DO
1359 0 : END SUBROUTINE allocate_dbcsr_matrix_set_3d
1360 :
1361 : ! **************************************************************************************************
1362 : !> \brief Allocate and initialize a real matrix 4-dimensional set.
1363 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1364 : !> \param[in] nmatrix Size of set
1365 : !> \param mmatrix ...
1366 : !> \param pmatrix ...
1367 : !> \par History
1368 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1369 : ! **************************************************************************************************
1370 0 : SUBROUTINE allocate_dbcsr_matrix_set_4d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix)
1371 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1372 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix
1373 :
1374 : INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1375 :
1376 0 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1377 0 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix))
1378 0 : DO lmatrix = 1, qmatrix
1379 0 : DO kmatrix = 1, pmatrix
1380 0 : DO jmatrix = 1, mmatrix
1381 0 : DO imatrix = 1, nmatrix
1382 0 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1383 : END DO
1384 : END DO
1385 : END DO
1386 : END DO
1387 0 : END SUBROUTINE allocate_dbcsr_matrix_set_4d
1388 :
1389 : ! **************************************************************************************************
1390 : !> \brief Allocate and initialize a real matrix 5-dimensional set.
1391 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1392 : !> \param[in] nmatrix Size of set
1393 : !> \param mmatrix ...
1394 : !> \param pmatrix ...
1395 : !> \par History
1396 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1397 : ! **************************************************************************************************
1398 0 : SUBROUTINE allocate_dbcsr_matrix_set_5d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix, smatrix)
1399 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), POINTER :: matrix_set
1400 : INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix, smatrix
1401 :
1402 : INTEGER :: imatrix, jmatrix, kmatrix, lmatrix, hmatrix
1403 :
1404 0 : IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1405 0 : ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix, smatrix))
1406 0 : DO hmatrix = 1, smatrix
1407 0 : DO lmatrix = 1, qmatrix
1408 0 : DO kmatrix = 1, pmatrix
1409 0 : DO jmatrix = 1, mmatrix
1410 0 : DO imatrix = 1, nmatrix
1411 0 : NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1412 : END DO
1413 : END DO
1414 : END DO
1415 : END DO
1416 : END DO
1417 0 : END SUBROUTINE allocate_dbcsr_matrix_set_5d
1418 :
1419 : ! **************************************************************************************************
1420 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1421 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1422 : !> \par History
1423 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1424 : ! **************************************************************************************************
1425 178728 : SUBROUTINE deallocate_dbcsr_matrix_set_1d(matrix_set)
1426 :
1427 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1428 :
1429 : INTEGER :: imatrix
1430 :
1431 178728 : IF (ASSOCIATED(matrix_set)) THEN
1432 479656 : DO imatrix = 1, SIZE(matrix_set)
1433 479656 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix)%matrix)
1434 : END DO
1435 177246 : DEALLOCATE (matrix_set)
1436 : END IF
1437 :
1438 178728 : END SUBROUTINE deallocate_dbcsr_matrix_set_1d
1439 :
1440 : ! **************************************************************************************************
1441 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1442 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1443 : !> \par History
1444 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1445 : ! **************************************************************************************************
1446 139185 : SUBROUTINE deallocate_dbcsr_matrix_set_2d(matrix_set)
1447 :
1448 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1449 :
1450 : INTEGER :: imatrix, jmatrix
1451 :
1452 139185 : IF (ASSOCIATED(matrix_set)) THEN
1453 741124 : DO jmatrix = 1, SIZE(matrix_set, 2)
1454 1453107 : DO imatrix = 1, SIZE(matrix_set, 1)
1455 1316924 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix)%matrix)
1456 : END DO
1457 : END DO
1458 136183 : DEALLOCATE (matrix_set)
1459 : END IF
1460 139185 : END SUBROUTINE deallocate_dbcsr_matrix_set_2d
1461 :
1462 : ! **************************************************************************************************
1463 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1464 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1465 : !> \par History
1466 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1467 : ! **************************************************************************************************
1468 0 : SUBROUTINE deallocate_dbcsr_matrix_set_3d(matrix_set)
1469 :
1470 : TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1471 :
1472 : INTEGER :: imatrix, jmatrix, kmatrix
1473 :
1474 0 : IF (ASSOCIATED(matrix_set)) THEN
1475 0 : DO kmatrix = 1, SIZE(matrix_set, 3)
1476 0 : DO jmatrix = 1, SIZE(matrix_set, 2)
1477 0 : DO imatrix = 1, SIZE(matrix_set, 1)
1478 0 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1479 : END DO
1480 : END DO
1481 : END DO
1482 0 : DEALLOCATE (matrix_set)
1483 : END IF
1484 0 : END SUBROUTINE deallocate_dbcsr_matrix_set_3d
1485 :
1486 : ! **************************************************************************************************
1487 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1488 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1489 : !> \par History
1490 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1491 : ! **************************************************************************************************
1492 0 : SUBROUTINE deallocate_dbcsr_matrix_set_4d(matrix_set)
1493 :
1494 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1495 :
1496 : INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1497 :
1498 0 : IF (ASSOCIATED(matrix_set)) THEN
1499 0 : DO lmatrix = 1, SIZE(matrix_set, 4)
1500 0 : DO kmatrix = 1, SIZE(matrix_set, 3)
1501 0 : DO jmatrix = 1, SIZE(matrix_set, 2)
1502 0 : DO imatrix = 1, SIZE(matrix_set, 1)
1503 0 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1504 : END DO
1505 : END DO
1506 : END DO
1507 : END DO
1508 0 : DEALLOCATE (matrix_set)
1509 : END IF
1510 0 : END SUBROUTINE deallocate_dbcsr_matrix_set_4d
1511 :
1512 : ! **************************************************************************************************
1513 : !> \brief Deallocate a real matrix set and release all of the member matrices.
1514 : !> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1515 : !> \par History
1516 : !> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1517 : ! **************************************************************************************************
1518 0 : SUBROUTINE deallocate_dbcsr_matrix_set_5d(matrix_set)
1519 :
1520 : TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), POINTER :: matrix_set
1521 :
1522 : INTEGER :: imatrix, jmatrix, kmatrix, hmatrix, lmatrix
1523 :
1524 0 : IF (ASSOCIATED(matrix_set)) THEN
1525 0 : DO hmatrix = 1, SIZE(matrix_set, 5)
1526 0 : DO lmatrix = 1, SIZE(matrix_set, 4)
1527 0 : DO kmatrix = 1, SIZE(matrix_set, 3)
1528 0 : DO jmatrix = 1, SIZE(matrix_set, 2)
1529 0 : DO imatrix = 1, SIZE(matrix_set, 1)
1530 0 : CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1531 : END DO
1532 : END DO
1533 : END DO
1534 : END DO
1535 : END DO
1536 0 : DEALLOCATE (matrix_set)
1537 : END IF
1538 0 : END SUBROUTINE deallocate_dbcsr_matrix_set_5d
1539 :
1540 : END MODULE cp_dbcsr_operations
|
