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 : MODULE qs_fb_filter_matrix_methods
9 :
10 : USE atomic_kind_types, ONLY: atomic_kind_type,&
11 : get_atomic_kind
12 : USE cp_dbcsr_api, ONLY: dbcsr_create,&
13 : dbcsr_distribution_type,&
14 : dbcsr_finalize,&
15 : dbcsr_get_info,&
16 : dbcsr_get_stored_coordinates,&
17 : dbcsr_put_block,&
18 : dbcsr_type,&
19 : dbcsr_type_no_symmetry
20 : USE input_constants, ONLY: smear_fermi_dirac
21 : USE kinds, ONLY: default_string_length,&
22 : dp,&
23 : int_8
24 : USE message_passing, ONLY: mp_para_env_type
25 : USE particle_types, ONLY: particle_type
26 : USE qs_fb_atomic_halo_types, ONLY: fb_atomic_halo_create,&
27 : fb_atomic_halo_get,&
28 : fb_atomic_halo_list_get,&
29 : fb_atomic_halo_list_obj,&
30 : fb_atomic_halo_nullify,&
31 : fb_atomic_halo_obj,&
32 : fb_atomic_halo_release,&
33 : fb_atomic_halo_set
34 : USE qs_fb_atomic_matrix_methods, ONLY: fb_atmatrix_calc_size,&
35 : fb_atmatrix_construct,&
36 : fb_atmatrix_construct_2,&
37 : fb_atmatrix_generate_com_pairs_2
38 : USE qs_fb_com_tasks_types, ONLY: &
39 : TASK_COST, TASK_DEST, TASK_N_RECORDS, TASK_PAIR, TASK_SRC, &
40 : fb_com_atom_pairs_calc_buffer_sizes, fb_com_atom_pairs_create, fb_com_atom_pairs_decode, &
41 : fb_com_atom_pairs_distribute_blks, fb_com_atom_pairs_gather_blks, fb_com_atom_pairs_get, &
42 : fb_com_atom_pairs_has_data, fb_com_atom_pairs_init, fb_com_atom_pairs_nullify, &
43 : fb_com_atom_pairs_obj, fb_com_atom_pairs_release, fb_com_tasks_build_atom_pairs, &
44 : fb_com_tasks_create, fb_com_tasks_encode_pair, fb_com_tasks_nullify, fb_com_tasks_obj, &
45 : fb_com_tasks_release, fb_com_tasks_set, fb_com_tasks_transpose_dest_src
46 : USE qs_fb_matrix_data_types, ONLY: fb_matrix_data_add,&
47 : fb_matrix_data_create,&
48 : fb_matrix_data_has_data,&
49 : fb_matrix_data_nullify,&
50 : fb_matrix_data_obj,&
51 : fb_matrix_data_release
52 : USE qs_fb_trial_fns_types, ONLY: fb_trial_fns_get,&
53 : fb_trial_fns_obj
54 : USE smearing_utils, ONLY: SmearFixed,&
55 : SmearOcc
56 : USE string_utilities, ONLY: compress,&
57 : uppercase
58 : #include "./base/base_uses.f90"
59 :
60 : IMPLICIT NONE
61 :
62 : PRIVATE
63 :
64 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_fb_filter_matrix_methods'
65 :
66 : PUBLIC :: fb_fltrmat_build, &
67 : fb_fltrmat_build_2
68 :
69 : CONTAINS
70 :
71 : ! **************************************************************************************************
72 : !> \brief Build the filter matrix, with MPI communications happening at each
73 : !> step. Less efficient on communication, but more efficient on
74 : !> memory usage (compared to fb_fltrmat_build_2)
75 : !> \param H_mat : DBCSR system KS matrix
76 : !> \param S_mat : DBCSR system overlap matrix
77 : !> \param atomic_halos : list of all local atomic halos, each halo gives
78 : !> one atomic matrix and contributes to one blk
79 : !> col to the filter matrix
80 : !> \param trial_fns : the trial functions to be used to shrink the
81 : !> size of the new "filtered" basis
82 : !> \param para_env : cp2k parallel environment
83 : !> \param particle_set : set of all particles in the system
84 : !> \param fermi_level : the fermi level used for defining the filter
85 : !> function, which is a Fermi-Dirac distribution
86 : !> function
87 : !> \param filter_temp : the filter temperature used for defining the
88 : !> filter function
89 : !> \param name : name given to the filter matrix
90 : !> \param filter_mat : DBCSR format filter matrix
91 : !> \param tolerance : anything less than tolerance is treated as zero
92 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
93 : ! **************************************************************************************************
94 64 : SUBROUTINE fb_fltrmat_build(H_mat, &
95 : S_mat, &
96 : atomic_halos, &
97 : trial_fns, &
98 : para_env, &
99 : particle_set, &
100 : fermi_level, &
101 : filter_temp, &
102 : name, &
103 : filter_mat, &
104 : tolerance)
105 : TYPE(dbcsr_type), POINTER :: H_mat, S_mat
106 : TYPE(fb_atomic_halo_list_obj), INTENT(IN) :: atomic_halos
107 : TYPE(fb_trial_fns_obj), INTENT(IN) :: trial_fns
108 : TYPE(mp_para_env_type), POINTER :: para_env
109 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
110 : REAL(KIND=dp), INTENT(IN) :: fermi_level, filter_temp
111 : CHARACTER(LEN=*), INTENT(IN) :: name
112 : TYPE(dbcsr_type), POINTER :: filter_mat
113 : REAL(KIND=dp), INTENT(IN) :: tolerance
114 :
115 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_build'
116 :
117 : CHARACTER(LEN=32) :: symmetry_string
118 : CHARACTER(LEN=default_string_length) :: name_string
119 : INTEGER :: handle, iblkcol, ihalo, ikind, &
120 : max_nhalos, nblkcols_total, nhalos
121 64 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, dummy_halo_atoms, ntfns, &
122 64 : row_blk_size
123 : LOGICAL :: send_data_only
124 : TYPE(atomic_kind_type), POINTER :: atomic_kind
125 : TYPE(dbcsr_distribution_type) :: dbcsr_dist
126 : TYPE(fb_atomic_halo_obj) :: dummy_atomic_halo
127 64 : TYPE(fb_atomic_halo_obj), DIMENSION(:), POINTER :: halos
128 :
129 64 : CALL timeset(routineN, handle)
130 :
131 64 : NULLIFY (halos, atomic_kind, ntfns, dummy_halo_atoms, row_blk_size, col_blk_size)
132 64 : CALL fb_atomic_halo_nullify(dummy_atomic_halo)
133 :
134 : ! filter_mat must be of a dissassociated status (i.e. brand new)
135 64 : CPASSERT(.NOT. ASSOCIATED(filter_mat))
136 :
137 : ! get trial function information
138 : CALL fb_trial_fns_get(trial_fns=trial_fns, &
139 64 : nfunctions=ntfns)
140 :
141 : ! calculate the row_blk_size and col_blk_size arrays for
142 : ! constructing the filter matrix in DBCSR format
143 : ! row_blk_size for the filter matrix is the same as H or S
144 : CALL dbcsr_get_info(H_mat, &
145 : nblkcols_total=nblkcols_total, &
146 : row_blk_size=row_blk_size, &
147 64 : distribution=dbcsr_dist)
148 192 : ALLOCATE (col_blk_size(nblkcols_total))
149 576 : col_blk_size = 0
150 576 : DO iblkcol = 1, nblkcols_total
151 512 : atomic_kind => particle_set(iblkcol)%atomic_kind
152 : CALL get_atomic_kind(atomic_kind=atomic_kind, &
153 512 : kind_number=ikind)
154 576 : col_blk_size(iblkcol) = ntfns(ikind)
155 : END DO
156 : ! DO NOT deallocate cbs if gift=.TRUE. as col_blk_sizes will only point to cbs
157 64 : name_string = name
158 64 : CALL compress(name_string)
159 64 : CALL uppercase(name_string)
160 : ! the filter matrix is non-square and is always non-symmetric
161 64 : symmetry_string = dbcsr_type_no_symmetry
162 : ! create empty filter matrix
163 64 : ALLOCATE (filter_mat)
164 : CALL dbcsr_create(matrix=filter_mat, &
165 : name=name_string, &
166 : dist=dbcsr_dist, &
167 : matrix_type=symmetry_string, &
168 : row_blk_size=row_blk_size, &
169 64 : col_blk_size=col_blk_size)
170 64 : DEALLOCATE (col_blk_size)
171 :
172 : CALL fb_atomic_halo_list_get(atomic_halos=atomic_halos, &
173 : nhalos=nhalos, &
174 : max_nhalos=max_nhalos, &
175 64 : halos=halos)
176 :
177 : ! create dummy empty atomic halo
178 64 : CALL fb_atomic_halo_create(dummy_atomic_halo)
179 64 : ALLOCATE (dummy_halo_atoms(0))
180 : CALL fb_atomic_halo_set(atomic_halo=dummy_atomic_halo, &
181 : owner_atom=0, &
182 : owner_id_in_halo=0, &
183 : natoms=0, &
184 : halo_atoms=dummy_halo_atoms, &
185 : nelectrons=0, &
186 64 : sorted=.TRUE.)
187 :
188 64 : send_data_only = .FALSE.
189 :
190 320 : DO ihalo = 1, max_nhalos
191 256 : IF (ihalo > nhalos) THEN
192 : send_data_only = .TRUE.
193 : END IF
194 : ! construct the filter matrix block by block
195 320 : IF (send_data_only) THEN
196 : CALL fb_fltrmat_add_blkcol(H_mat, &
197 : S_mat, &
198 : dummy_atomic_halo, &
199 : trial_fns, &
200 : para_env, &
201 : particle_set, &
202 : fermi_level, &
203 : filter_temp, &
204 : filter_mat, &
205 0 : tolerance)
206 : ELSE
207 : CALL fb_fltrmat_add_blkcol(H_mat, &
208 : S_mat, &
209 : halos(ihalo), &
210 : trial_fns, &
211 : para_env, &
212 : particle_set, &
213 : fermi_level, &
214 : filter_temp, &
215 : filter_mat, &
216 256 : tolerance)
217 : END IF ! send_data_only
218 : END DO
219 :
220 : ! finalise the filter matrix
221 64 : CALL dbcsr_finalize(filter_mat)
222 :
223 : ! cleanup
224 64 : CALL fb_atomic_halo_release(dummy_atomic_halo)
225 :
226 64 : CALL timestop(handle)
227 :
228 192 : END SUBROUTINE fb_fltrmat_build
229 :
230 : ! **************************************************************************************************
231 : !> \brief Build the filter matrix, with MPI communications grouped together.
232 : !> More effcient on communication, less efficient on memory (compared
233 : !> to fb_fltrmat_build)
234 : !> \param H_mat : DBCSR system KS matrix
235 : !> \param S_mat : DBCSR system overlap matrix
236 : !> \param atomic_halos : list of all local atomic halos, each halo gives
237 : !> one atomic matrix and contributes to one blk
238 : !> col to the filter matrix
239 : !> \param trial_fns : the trial functions to be used to shrink the
240 : !> size of the new "filtered" basis
241 : !> \param para_env : cp2k parallel environment
242 : !> \param particle_set : set of all particles in the system
243 : !> \param fermi_level : the fermi level used for defining the filter
244 : !> function, which is a Fermi-Dirac distribution
245 : !> function
246 : !> \param filter_temp : the filter temperature used for defining the
247 : !> filter function
248 : !> \param name : name given to the filter matrix
249 : !> \param filter_mat : DBCSR format filter matrix
250 : !> \param tolerance : anything less than tolerance is treated as zero
251 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
252 : ! **************************************************************************************************
253 16 : SUBROUTINE fb_fltrmat_build_2(H_mat, &
254 : S_mat, &
255 : atomic_halos, &
256 : trial_fns, &
257 : para_env, &
258 : particle_set, &
259 : fermi_level, &
260 : filter_temp, &
261 : name, &
262 : filter_mat, &
263 : tolerance)
264 : TYPE(dbcsr_type), POINTER :: H_mat, S_mat
265 : TYPE(fb_atomic_halo_list_obj), INTENT(IN) :: atomic_halos
266 : TYPE(fb_trial_fns_obj), INTENT(IN) :: trial_fns
267 : TYPE(mp_para_env_type), POINTER :: para_env
268 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
269 : REAL(KIND=dp), INTENT(IN) :: fermi_level, filter_temp
270 : CHARACTER(LEN=*), INTENT(IN) :: name
271 : TYPE(dbcsr_type), POINTER :: filter_mat
272 : REAL(KIND=dp), INTENT(IN) :: tolerance
273 :
274 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_build_2'
275 :
276 : CHARACTER(LEN=default_string_length) :: name_string
277 : INTEGER :: handle, iblkcol, ihalo, ikind, &
278 : natoms_global, natoms_in_halo, &
279 : nblkcols_total, nblks_recv, nhalos, &
280 : nmax
281 16 : INTEGER, DIMENSION(:), POINTER :: col_blk_size, ntfns, row_blk_size
282 : LOGICAL :: check_ok
283 : TYPE(atomic_kind_type), POINTER :: atomic_kind
284 : TYPE(dbcsr_distribution_type) :: dbcsr_dist
285 16 : TYPE(fb_atomic_halo_obj), DIMENSION(:), POINTER :: halos
286 : TYPE(fb_com_atom_pairs_obj) :: atmatrix_blks_recv, atmatrix_blks_send, &
287 : filter_mat_blks_recv, &
288 : filter_mat_blks_send
289 : TYPE(fb_matrix_data_obj) :: filter_mat_data, H_mat_data, S_mat_data
290 :
291 16 : CALL timeset(routineN, handle)
292 :
293 16 : NULLIFY (halos, atomic_kind, row_blk_size, col_blk_size, ntfns)
294 :
295 : ! filter_mat must be of a dissassociated status (i.e. brand new)
296 16 : check_ok = .NOT. ASSOCIATED(filter_mat)
297 16 : CPASSERT(check_ok)
298 :
299 : ! get total number of atoms
300 16 : natoms_global = SIZE(particle_set)
301 :
302 : ! get trial function information
303 : CALL fb_trial_fns_get(trial_fns=trial_fns, &
304 16 : nfunctions=ntfns)
305 :
306 : ! calculate the row_blk_size and col_blk_size arrays for
307 : ! constructing the filter matrix in DBCSR format
308 : ! row_blk_size for the filter matrix is the same as H or S
309 : CALL dbcsr_get_info(H_mat, &
310 : nblkcols_total=nblkcols_total, &
311 : row_blk_size=row_blk_size, &
312 16 : distribution=dbcsr_dist)
313 48 : ALLOCATE (col_blk_size(nblkcols_total))
314 144 : col_blk_size = 0
315 144 : DO iblkcol = 1, nblkcols_total
316 128 : atomic_kind => particle_set(iblkcol)%atomic_kind
317 : CALL get_atomic_kind(atomic_kind=atomic_kind, &
318 128 : kind_number=ikind)
319 144 : col_blk_size(iblkcol) = ntfns(ikind)
320 : END DO
321 : ! DO NOT deallocate cbs if gift=.TRUE. as col_blk_sizes will only point to cbs
322 16 : name_string = name
323 16 : CALL compress(name_string)
324 16 : CALL uppercase(name_string)
325 : ! create empty filter matrix (it is always non-symmetric as it is non-square)
326 16 : ALLOCATE (filter_mat)
327 : CALL dbcsr_create(matrix=filter_mat, &
328 : name=name_string, &
329 : dist=dbcsr_dist, &
330 : matrix_type=dbcsr_type_no_symmetry, &
331 : row_blk_size=row_blk_size, &
332 16 : col_blk_size=col_blk_size)
333 16 : DEALLOCATE (col_blk_size)
334 :
335 : ! get all the blocks required for constructing atomic matrics, and
336 : ! store it in a fb_matrix_data object
337 16 : CALL fb_matrix_data_nullify(H_mat_data)
338 16 : CALL fb_matrix_data_nullify(S_mat_data)
339 16 : CALL fb_com_atom_pairs_nullify(atmatrix_blks_send)
340 16 : CALL fb_com_atom_pairs_nullify(atmatrix_blks_recv)
341 16 : CALL fb_com_atom_pairs_create(atmatrix_blks_send)
342 16 : CALL fb_com_atom_pairs_create(atmatrix_blks_recv)
343 : ! H matrix
344 : CALL fb_atmatrix_generate_com_pairs_2(H_mat, &
345 : atomic_halos, &
346 : para_env, &
347 : atmatrix_blks_send, &
348 16 : atmatrix_blks_recv)
349 : CALL fb_com_atom_pairs_get(atom_pairs=atmatrix_blks_recv, &
350 16 : npairs=nblks_recv)
351 : CALL fb_matrix_data_create(H_mat_data, &
352 : nblks_recv, &
353 16 : natoms_global)
354 : CALL fb_com_atom_pairs_gather_blks(H_mat, &
355 : atmatrix_blks_send, &
356 : atmatrix_blks_recv, &
357 : para_env, &
358 16 : H_mat_data)
359 : ! S matrix
360 : CALL fb_atmatrix_generate_com_pairs_2(S_mat, &
361 : atomic_halos, &
362 : para_env, &
363 : atmatrix_blks_send, &
364 16 : atmatrix_blks_recv)
365 : CALL fb_com_atom_pairs_get(atom_pairs=atmatrix_blks_recv, &
366 16 : npairs=nblks_recv)
367 : CALL fb_matrix_data_create(S_mat_data, &
368 : nblks_recv, &
369 16 : natoms_global)
370 : CALL fb_com_atom_pairs_gather_blks(S_mat, &
371 : atmatrix_blks_send, &
372 : atmatrix_blks_recv, &
373 : para_env, &
374 16 : S_mat_data)
375 : ! cleanup
376 16 : CALL fb_com_atom_pairs_release(atmatrix_blks_send)
377 16 : CALL fb_com_atom_pairs_release(atmatrix_blks_recv)
378 :
379 : ! make filter matrix blocks one by one and store in an
380 : ! matrix_data_obj
381 16 : CALL fb_matrix_data_nullify(filter_mat_data)
382 : CALL fb_atomic_halo_list_get(atomic_halos=atomic_halos, &
383 : nhalos=nhalos, &
384 16 : halos=halos)
385 16 : nmax = 0
386 80 : DO ihalo = 1, nhalos
387 : CALL fb_atomic_halo_get(atomic_halo=halos(ihalo), &
388 64 : natoms=natoms_in_halo)
389 80 : nmax = nmax + natoms_in_halo
390 : END DO
391 : CALL fb_matrix_data_create(filter_mat_data, &
392 : nmax, &
393 16 : natoms_global)
394 80 : DO ihalo = 1, nhalos
395 : CALL fb_fltrmat_add_blkcol_2(H_mat, &
396 : S_mat, &
397 : H_mat_data, &
398 : S_mat_data, &
399 : halos(ihalo), &
400 : trial_fns, &
401 : particle_set, &
402 : fermi_level, &
403 : filter_temp, &
404 : filter_mat_data, &
405 80 : tolerance)
406 : END DO
407 : ! clean up
408 16 : CALL fb_matrix_data_release(H_mat_data)
409 16 : CALL fb_matrix_data_release(S_mat_data)
410 :
411 : ! distribute the relevant blocks from the matrix_data_obj to DBCSR
412 : ! filter matrix
413 16 : CALL fb_com_atom_pairs_nullify(filter_mat_blks_send)
414 16 : CALL fb_com_atom_pairs_nullify(filter_mat_blks_recv)
415 16 : CALL fb_com_atom_pairs_create(filter_mat_blks_send)
416 16 : CALL fb_com_atom_pairs_create(filter_mat_blks_recv)
417 : CALL fb_fltrmat_generate_com_pairs_2(filter_mat, &
418 : atomic_halos, &
419 : para_env, &
420 : filter_mat_blks_send, &
421 16 : filter_mat_blks_recv)
422 : CALL fb_com_atom_pairs_distribute_blks(filter_mat_data, &
423 : filter_mat_blks_send, &
424 : filter_mat_blks_recv, &
425 : para_env, &
426 16 : filter_mat)
427 : ! cleanup
428 16 : CALL fb_com_atom_pairs_release(filter_mat_blks_send)
429 16 : CALL fb_com_atom_pairs_release(filter_mat_blks_recv)
430 16 : CALL fb_matrix_data_release(filter_mat_data)
431 :
432 : ! finalise matrix
433 16 : CALL dbcsr_finalize(filter_mat)
434 :
435 16 : CALL timestop(handle)
436 :
437 96 : END SUBROUTINE fb_fltrmat_build_2
438 :
439 : ! **************************************************************************************************
440 : !> \brief Add a computed blocks in one column to the filter matrix. This
441 : !> version is used by fb_fltrmat_build, for the case where MPI
442 : !> communications are done at each step
443 : !> It does not finalise the filter matrix
444 : !> \param H_mat : DBCSR system KS matrix
445 : !> \param S_mat : DBCSR system overlap matrix
446 : !> \param atomic_halo : the halo that contributes to the blk
447 : !> col of the filter matrix
448 : !> \param trial_fns ...
449 : !> \param para_env : cp2k parallel environment
450 : !> \param particle_set : set of all particles in the system
451 : !> \param fermi_level : the fermi level used for defining the filter
452 : !> function, which is a Fermi-Dirac distribution
453 : !> function
454 : !> \param filter_temp : the filter temperature used for defining the
455 : !> filter function
456 : !> \param filter_mat : DBCSR format filter matrix
457 : !> \param tolerance : anything smaller than tolerance is treated as zero
458 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
459 : ! **************************************************************************************************
460 256 : SUBROUTINE fb_fltrmat_add_blkcol(H_mat, &
461 : S_mat, &
462 : atomic_halo, &
463 : trial_fns, &
464 : para_env, &
465 : particle_set, &
466 : fermi_level, &
467 : filter_temp, &
468 : filter_mat, &
469 : tolerance)
470 : TYPE(dbcsr_type), POINTER :: H_mat, S_mat
471 : TYPE(fb_atomic_halo_obj), INTENT(IN) :: atomic_halo
472 : TYPE(fb_trial_fns_obj), INTENT(IN) :: trial_fns
473 : TYPE(mp_para_env_type), POINTER :: para_env
474 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
475 : REAL(KIND=dp), INTENT(IN) :: fermi_level, filter_temp
476 : TYPE(dbcsr_type), POINTER :: filter_mat
477 : REAL(KIND=dp), INTENT(IN) :: tolerance
478 :
479 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_add_blkcol'
480 :
481 : INTEGER :: handle, handle_mpi, iatom_global, iatom_in_halo, ind, ipair, ipe, itrial, &
482 : jatom_global, jatom_in_halo, jkind, natoms_global, natoms_in_halo, ncols_atmatrix, &
483 : ncols_blk, nrows_atmatrix, nrows_blk, numprocs, pe, recv_encode, send_encode
484 256 : INTEGER(KIND=int_8), DIMENSION(:), POINTER :: pairs_recv, pairs_send
485 256 : INTEGER, ALLOCATABLE, DIMENSION(:) :: atomic_H_blk_col_start, atomic_H_blk_row_start, &
486 256 : atomic_S_blk_col_start, atomic_S_blk_row_start, col_block_size_data, ind_in_halo, &
487 256 : recv_disps, recv_pair_count, recv_pair_disps, recv_sizes, send_disps, send_pair_count, &
488 256 : send_pair_disps, send_sizes
489 256 : INTEGER, DIMENSION(:), POINTER :: halo_atoms, ntfns, row_block_size_data
490 256 : INTEGER, DIMENSION(:, :), POINTER :: tfns
491 256 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), TARGET :: recv_buf, send_buf
492 256 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: atomic_filter_mat, atomic_H, atomic_S
493 256 : REAL(kind=dp), DIMENSION(:), POINTER :: vector
494 : TYPE(atomic_kind_type), POINTER :: atomic_kind
495 : TYPE(fb_com_atom_pairs_obj) :: com_pairs_recv, com_pairs_send
496 :
497 256 : CALL timeset(routineN, handle)
498 :
499 256 : NULLIFY (atomic_kind, halo_atoms, ntfns, pairs_send, pairs_recv, &
500 256 : row_block_size_data, tfns)
501 256 : CALL fb_com_atom_pairs_nullify(com_pairs_send)
502 256 : CALL fb_com_atom_pairs_nullify(com_pairs_recv)
503 :
504 : ! ----------------------------------------------------------------------
505 : ! Get communication buffers ready
506 : ! ----------------------------------------------------------------------
507 :
508 : ! generate send and recv atom pairs
509 256 : CALL fb_com_atom_pairs_create(com_pairs_send)
510 256 : CALL fb_com_atom_pairs_create(com_pairs_recv)
511 : CALL fb_fltrmat_generate_com_pairs(filter_mat, &
512 : atomic_halo, &
513 : para_env, &
514 : com_pairs_send, &
515 256 : com_pairs_recv)
516 : CALL fb_com_atom_pairs_get(atom_pairs=com_pairs_send, &
517 : natoms_encode=send_encode, &
518 256 : pairs=pairs_send)
519 : CALL fb_com_atom_pairs_get(atom_pairs=com_pairs_recv, &
520 : natoms_encode=recv_encode, &
521 256 : pairs=pairs_recv)
522 :
523 : ! get para_env info
524 256 : numprocs = para_env%num_pe
525 : ! me = para_env%mepos + 1 ! my process id, starting counting from 1
526 :
527 : ! obtain trail function information
528 : CALL fb_trial_fns_get(trial_fns=trial_fns, &
529 : nfunctions=ntfns, &
530 256 : functions=tfns)
531 :
532 : ! obtain row and col block size data for filter matrix
533 256 : CALL dbcsr_get_info(H_mat, row_blk_size=row_block_size_data)
534 256 : natoms_global = SIZE(particle_set)
535 768 : ALLOCATE (col_block_size_data(natoms_global))
536 2304 : DO jatom_global = 1, natoms_global
537 2048 : atomic_kind => particle_set(jatom_global)%atomic_kind
538 2048 : CALL get_atomic_kind(atomic_kind=atomic_kind, kind_number=jkind)
539 2304 : col_block_size_data(jatom_global) = ntfns(jkind)
540 : END DO
541 :
542 : ! allocate temporary arrays for send
543 768 : ALLOCATE (send_sizes(numprocs))
544 512 : ALLOCATE (send_disps(numprocs))
545 512 : ALLOCATE (send_pair_count(numprocs))
546 512 : ALLOCATE (send_pair_disps(numprocs))
547 : ! setup send buffer sizes
548 : CALL fb_com_atom_pairs_calc_buffer_sizes(com_pairs_send, &
549 : numprocs, &
550 : row_block_size_data, &
551 : col_block_size_data, &
552 : send_sizes, &
553 : send_disps, &
554 : send_pair_count, &
555 256 : send_pair_disps)
556 : ! allocate send buffer
557 1280 : ALLOCATE (send_buf(SUM(send_sizes)))
558 :
559 : ! allocate temporary array for recv
560 512 : ALLOCATE (recv_sizes(numprocs))
561 512 : ALLOCATE (recv_disps(numprocs))
562 512 : ALLOCATE (recv_pair_count(numprocs))
563 512 : ALLOCATE (recv_pair_disps(numprocs))
564 : ! setup recv buffer sizes
565 : CALL fb_com_atom_pairs_calc_buffer_sizes(com_pairs_recv, &
566 : numprocs, &
567 : row_block_size_data, &
568 : col_block_size_data, &
569 : recv_sizes, &
570 : recv_disps, &
571 : recv_pair_count, &
572 256 : recv_pair_disps)
573 : ! allocate recv buffer
574 1280 : ALLOCATE (recv_buf(SUM(recv_sizes)))
575 :
576 : ! ----------------------------------------------------------------------
577 : ! Construct atomic filter matrix for this atomic_halo
578 : ! ----------------------------------------------------------------------
579 :
580 : CALL fb_atomic_halo_get(atomic_halo=atomic_halo, &
581 : natoms=natoms_in_halo, &
582 256 : halo_atoms=halo_atoms)
583 :
584 : ! construct atomic matrix for H for atomic_halo
585 : ALLOCATE (atomic_H_blk_row_start(natoms_in_halo + 1), &
586 1024 : atomic_H_blk_col_start(natoms_in_halo + 1))
587 : CALL fb_atmatrix_calc_size(H_mat, &
588 : atomic_halo, &
589 : nrows_atmatrix, &
590 : ncols_atmatrix, &
591 : atomic_H_blk_row_start, &
592 256 : atomic_H_blk_col_start)
593 :
594 1024 : ALLOCATE (atomic_H(nrows_atmatrix, ncols_atmatrix))
595 : CALL fb_atmatrix_construct(H_mat, &
596 : atomic_halo, &
597 : para_env, &
598 : atomic_H, &
599 : atomic_H_blk_row_start, &
600 256 : atomic_H_blk_col_start)
601 :
602 : ! construct atomic matrix for S for atomic_halo
603 : ALLOCATE (atomic_S_blk_row_start(natoms_in_halo + 1), &
604 768 : atomic_S_blk_col_start(natoms_in_halo + 1))
605 : CALL fb_atmatrix_calc_size(S_mat, &
606 : atomic_halo, &
607 : nrows_atmatrix, &
608 : ncols_atmatrix, &
609 : atomic_S_blk_row_start, &
610 256 : atomic_S_blk_col_start)
611 1024 : ALLOCATE (atomic_S(nrows_atmatrix, ncols_atmatrix))
612 : CALL fb_atmatrix_construct(S_mat, &
613 : atomic_halo, &
614 : para_env, &
615 : atomic_S, &
616 : atomic_S_blk_row_start, &
617 256 : atomic_S_blk_col_start)
618 :
619 : ! construct the atomic filter matrix
620 768 : ALLOCATE (atomic_filter_mat(nrows_atmatrix, ncols_atmatrix))
621 : ! calculate atomic filter matrix only if it is non-zero sized
622 256 : IF (nrows_atmatrix > 0 .AND. ncols_atmatrix > 0) THEN
623 : CALL fb_fltrmat_build_atomic_fltrmat(atomic_H, &
624 : atomic_S, &
625 : fermi_level, &
626 : filter_temp, &
627 : atomic_filter_mat, &
628 256 : tolerance)
629 : END IF
630 :
631 : ! ----------------------------------------------------------------------
632 : ! Construct filter matrix blocks and add to the correct locations
633 : ! in send_buffer
634 : ! ----------------------------------------------------------------------
635 :
636 : ! preconstruct iatom_global to iatom_in_halo map
637 512 : ALLOCATE (ind_in_halo(natoms_global))
638 256 : ind_in_halo = 0
639 2304 : DO iatom_in_halo = 1, natoms_in_halo
640 2048 : iatom_global = halo_atoms(iatom_in_halo)
641 2304 : ind_in_halo(iatom_global) = iatom_in_halo
642 : END DO
643 :
644 : ! initialise send buffer
645 256 : IF (SIZE(send_buf) > 0) send_buf = 0.0_dp
646 : ! assign values
647 768 : DO ipe = 1, numprocs
648 512 : send_sizes(ipe) = 0
649 2816 : DO ipair = 1, send_pair_count(ipe)
650 : CALL fb_com_atom_pairs_decode(pairs_send(send_pair_disps(ipe) + ipair), &
651 : pe, iatom_global, jatom_global, &
652 2048 : send_encode)
653 2048 : iatom_in_halo = ind_in_halo(iatom_global)
654 2048 : CPASSERT(iatom_in_halo > 0)
655 2048 : jatom_in_halo = ind_in_halo(jatom_global)
656 2048 : CPASSERT(jatom_in_halo > 0)
657 2048 : atomic_kind => particle_set(jatom_global)%atomic_kind
658 : CALL get_atomic_kind(atomic_kind=atomic_kind, &
659 2048 : kind_number=jkind)
660 2048 : nrows_blk = row_block_size_data(iatom_global)
661 2048 : ncols_blk = ntfns(jkind)
662 :
663 : ! do it column-wise one trial function at a time
664 10240 : DO itrial = 1, ntfns(jkind)
665 8192 : ind = send_disps(ipe) + send_sizes(ipe) + (itrial - 1)*nrows_blk
666 : CALL dgemv("N", &
667 : nrows_blk, &
668 : ncols_atmatrix, &
669 : 1.0_dp, &
670 : atomic_filter_mat( &
671 : atomic_H_blk_row_start(iatom_in_halo): &
672 : atomic_H_blk_row_start(iatom_in_halo + 1) - 1, &
673 : 1:ncols_atmatrix &
674 : ), &
675 : nrows_blk, &
676 : atomic_S( &
677 : 1:nrows_atmatrix, &
678 : atomic_S_blk_col_start(jatom_in_halo) + &
679 : tfns(itrial, jkind) - 1 &
680 : ), &
681 : 1, &
682 : 0.0_dp, &
683 : send_buf(ind + 1:ind + nrows_blk), &
684 23865344 : 1)
685 : END DO ! itrial
686 6656 : send_sizes(ipe) = send_sizes(ipe) + nrows_blk*ncols_blk
687 : END DO ! ipair
688 : END DO ! ipe
689 :
690 256 : DEALLOCATE (atomic_H)
691 256 : DEALLOCATE (atomic_H_blk_row_start)
692 256 : DEALLOCATE (atomic_S)
693 256 : DEALLOCATE (atomic_S_blk_row_start)
694 256 : DEALLOCATE (atomic_filter_mat)
695 256 : DEALLOCATE (ind_in_halo)
696 :
697 : ! ----------------------------------------------------------------------
698 : ! Do communication
699 : ! ----------------------------------------------------------------------
700 :
701 256 : CALL timeset("fb_fltrmat_add_blkcol_mpi", handle_mpi)
702 :
703 : CALL para_env%alltoall(send_buf, send_sizes, send_disps, &
704 256 : recv_buf, recv_sizes, recv_disps)
705 :
706 256 : CALL timestop(handle_mpi)
707 :
708 256 : DEALLOCATE (send_buf)
709 256 : DEALLOCATE (send_sizes)
710 256 : DEALLOCATE (send_disps)
711 256 : DEALLOCATE (send_pair_count)
712 256 : DEALLOCATE (send_pair_disps)
713 :
714 : ! ----------------------------------------------------------------------
715 : ! Unpack the recv buffer and add the blocks to correct parts of
716 : ! the DBCSR filter matrix
717 : ! ----------------------------------------------------------------------
718 :
719 768 : DO ipe = 1, numprocs
720 512 : recv_sizes(ipe) = 0
721 2816 : DO ipair = 1, recv_pair_count(ipe)
722 : CALL fb_com_atom_pairs_decode(pairs_recv(recv_pair_disps(ipe) + ipair), &
723 : pe, iatom_global, jatom_global, &
724 2048 : recv_encode)
725 2048 : nrows_blk = row_block_size_data(iatom_global)
726 2048 : ncols_blk = col_block_size_data(jatom_global)
727 2048 : ind = recv_disps(ipe) + recv_sizes(ipe)
728 2048 : vector => recv_buf((ind + 1):(ind + nrows_blk*ncols_blk))
729 : CALL dbcsr_put_block(filter_mat, &
730 : iatom_global, jatom_global, &
731 6144 : block=RESHAPE(vector, [nrows_blk, ncols_blk]))
732 4608 : recv_sizes(ipe) = recv_sizes(ipe) + nrows_blk*ncols_blk
733 : END DO ! ipair
734 : END DO ! ipe
735 :
736 : ! cleanup rest of the temporary arrays
737 256 : DEALLOCATE (recv_buf)
738 256 : DEALLOCATE (recv_sizes)
739 256 : DEALLOCATE (recv_pair_count)
740 256 : DEALLOCATE (recv_pair_disps)
741 :
742 256 : CALL fb_com_atom_pairs_release(com_pairs_send)
743 256 : CALL fb_com_atom_pairs_release(com_pairs_recv)
744 :
745 : ! cannot finalise the matrix until all blocks has been added
746 :
747 256 : CALL timestop(handle)
748 :
749 1792 : END SUBROUTINE fb_fltrmat_add_blkcol
750 :
751 : ! **************************************************************************************************
752 : !> \brief Computed blocks in one filter matrix column. This version is used by
753 : !> fb_fltrmat_build_2, where MPI communication is done collectively
754 : !> \param H_mat : DBCSR system KS matrix
755 : !> \param S_mat : DBCSR system overlap matrix
756 : !> \param H_mat_data : local storage of the relevant H_mat matrix blocks
757 : !> \param S_mat_data : local storage of the relevant S_mat matrix blocks
758 : !> \param atomic_halo : the halo that contributes to the blk
759 : !> col of the filter matrix
760 : !> \param trial_fns : trial functions data
761 : !> \param particle_set : set of all particles in the system
762 : !> \param fermi_level : the fermi level used for defining the filter
763 : !> function, which is a Fermi-Dirac distribution
764 : !> function
765 : !> \param filter_temp : the filter temperature used for defining the
766 : !> filter function
767 : !> \param filter_mat_data : local storage for the the computed filter matrix
768 : !> blocks
769 : !> \param tolerance : anything less than this is regarded as zero
770 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
771 : ! **************************************************************************************************
772 64 : SUBROUTINE fb_fltrmat_add_blkcol_2(H_mat, &
773 : S_mat, &
774 : H_mat_data, &
775 : S_mat_data, &
776 : atomic_halo, &
777 : trial_fns, &
778 : particle_set, &
779 : fermi_level, &
780 : filter_temp, &
781 : filter_mat_data, &
782 : tolerance)
783 : TYPE(dbcsr_type), POINTER :: H_mat, S_mat
784 : TYPE(fb_matrix_data_obj), INTENT(IN) :: H_mat_data, S_mat_data
785 : TYPE(fb_atomic_halo_obj), INTENT(IN) :: atomic_halo
786 : TYPE(fb_trial_fns_obj), INTENT(IN) :: trial_fns
787 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
788 : REAL(KIND=dp), INTENT(IN) :: fermi_level, filter_temp
789 : TYPE(fb_matrix_data_obj), INTENT(INOUT) :: filter_mat_data
790 : REAL(KIND=dp), INTENT(IN) :: tolerance
791 :
792 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_add_blkcol_2'
793 :
794 : INTEGER :: handle, iatom_global, iatom_in_halo, itrial, jatom_global, jatom_in_halo, jkind, &
795 : natoms_global, natoms_in_halo, ncols_atmatrix, ncols_blk, ncols_blk_max, nrows_atmatrix, &
796 : nrows_blk, nrows_blk_max
797 64 : INTEGER, ALLOCATABLE, DIMENSION(:) :: atomic_H_blk_col_start, atomic_H_blk_row_start, &
798 64 : atomic_S_blk_col_start, atomic_S_blk_row_start, col_block_size_data
799 64 : INTEGER, DIMENSION(:), POINTER :: halo_atoms, ntfns, row_block_size_data
800 64 : INTEGER, DIMENSION(:, :), POINTER :: tfns
801 : LOGICAL :: check_ok
802 64 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: atomic_filter_mat, atomic_H, atomic_S, &
803 64 : mat_blk
804 : TYPE(atomic_kind_type), POINTER :: atomic_kind
805 :
806 64 : CALL timeset(routineN, handle)
807 :
808 64 : NULLIFY (atomic_kind, halo_atoms, ntfns, row_block_size_data, tfns)
809 :
810 64 : check_ok = fb_matrix_data_has_data(H_mat_data)
811 64 : CPASSERT(check_ok)
812 64 : check_ok = fb_matrix_data_has_data(S_mat_data)
813 64 : CPASSERT(check_ok)
814 :
815 : ! obtain trial function information
816 : CALL fb_trial_fns_get(trial_fns=trial_fns, &
817 : nfunctions=ntfns, &
818 64 : functions=tfns)
819 :
820 : ! obtain row and col block size data for filter matrix
821 64 : CALL dbcsr_get_info(H_mat, row_blk_size=row_block_size_data)
822 64 : natoms_global = SIZE(particle_set)
823 192 : ALLOCATE (col_block_size_data(natoms_global))
824 576 : DO jatom_global = 1, natoms_global
825 512 : atomic_kind => particle_set(jatom_global)%atomic_kind
826 512 : CALL get_atomic_kind(atomic_kind=atomic_kind, kind_number=jkind)
827 576 : col_block_size_data(jatom_global) = ntfns(jkind)
828 : END DO
829 :
830 : ! ----------------------------------------------------------------------
831 : ! Construct atomic filter matrix for this atomic_halo
832 : ! ----------------------------------------------------------------------
833 :
834 : CALL fb_atomic_halo_get(atomic_halo=atomic_halo, &
835 : natoms=natoms_in_halo, &
836 64 : halo_atoms=halo_atoms)
837 :
838 : ! construct atomic matrix for H for atomic_halo
839 : ALLOCATE (atomic_H_blk_row_start(natoms_in_halo + 1), &
840 256 : atomic_H_blk_col_start(natoms_in_halo + 1))
841 : CALL fb_atmatrix_calc_size(H_mat, &
842 : atomic_halo, &
843 : nrows_atmatrix, &
844 : ncols_atmatrix, &
845 : atomic_H_blk_row_start, &
846 64 : atomic_H_blk_col_start)
847 256 : ALLOCATE (atomic_H(nrows_atmatrix, ncols_atmatrix))
848 : CALL fb_atmatrix_construct_2(H_mat_data, &
849 : atomic_halo, &
850 : atomic_H, &
851 : atomic_H_blk_row_start, &
852 64 : atomic_H_blk_col_start)
853 :
854 : ! construct atomic matrix for S for atomic_halo
855 : ALLOCATE (atomic_S_blk_row_start(natoms_in_halo + 1), &
856 192 : atomic_S_blk_col_start(natoms_in_halo + 1))
857 : CALL fb_atmatrix_calc_size(S_mat, &
858 : atomic_halo, &
859 : nrows_atmatrix, &
860 : ncols_atmatrix, &
861 : atomic_S_blk_row_start, &
862 64 : atomic_S_blk_col_start)
863 256 : ALLOCATE (atomic_S(nrows_atmatrix, ncols_atmatrix))
864 : CALL fb_atmatrix_construct_2(S_mat_data, &
865 : atomic_halo, &
866 : atomic_S, &
867 : atomic_S_blk_row_start, &
868 64 : atomic_S_blk_col_start)
869 :
870 : ! construct the atomic filter matrix
871 192 : ALLOCATE (atomic_filter_mat(nrows_atmatrix, ncols_atmatrix))
872 : ! calculate atomic filter matrix only if it is non-zero sized
873 64 : IF (nrows_atmatrix > 0 .AND. ncols_atmatrix > 0) THEN
874 : CALL fb_fltrmat_build_atomic_fltrmat(atomic_H, &
875 : atomic_S, &
876 : fermi_level, &
877 : filter_temp, &
878 : atomic_filter_mat, &
879 64 : tolerance)
880 : END IF
881 :
882 : ! ----------------------------------------------------------------------
883 : ! Construct filter matrix block and add to filter_mat_data
884 : ! ----------------------------------------------------------------------
885 :
886 : CALL fb_atomic_halo_get(atomic_halo=atomic_halo, &
887 : owner_atom=jatom_global, &
888 64 : owner_id_in_halo=jatom_in_halo)
889 576 : nrows_blk_max = MAXVAL(row_block_size_data)
890 128 : ncols_blk_max = MAXVAL(ntfns)
891 256 : ALLOCATE (mat_blk(nrows_blk_max, ncols_blk_max))
892 64 : mat_blk(:, :) = 0.0_dp
893 576 : DO iatom_in_halo = 1, natoms_in_halo
894 512 : iatom_global = halo_atoms(iatom_in_halo)
895 512 : atomic_kind => particle_set(jatom_global)%atomic_kind
896 : CALL get_atomic_kind(atomic_kind=atomic_kind, &
897 512 : kind_number=jkind)
898 512 : nrows_blk = row_block_size_data(iatom_global)
899 512 : ncols_blk = ntfns(jkind)
900 :
901 : ! do it column-wise one trial function at a time
902 2560 : DO itrial = 1, ntfns(jkind)
903 : CALL dgemv("N", &
904 : nrows_blk, &
905 : ncols_atmatrix, &
906 : 1.0_dp, &
907 : atomic_filter_mat( &
908 : atomic_H_blk_row_start(iatom_in_halo): &
909 : atomic_H_blk_row_start(iatom_in_halo + 1) - 1, &
910 : 1:ncols_atmatrix &
911 : ), &
912 : nrows_blk, &
913 : atomic_S( &
914 : 1:nrows_atmatrix, &
915 : atomic_S_blk_col_start(jatom_in_halo) + &
916 : tfns(itrial, jkind) - 1 &
917 : ), &
918 : 1, &
919 : 0.0_dp, &
920 : mat_blk( &
921 : 1:nrows_blk, &
922 : itrial), &
923 5966336 : 1)
924 : END DO ! itrial
925 : CALL fb_matrix_data_add(filter_mat_data, &
926 : iatom_global, &
927 : jatom_global, &
928 1088 : mat_blk(1:nrows_blk, 1:ncols_blk))
929 : ! DEALLOCATE(mat_blk, STAT=stat)
930 : ! CPPostcondition(stat==0, cp_failure_level, routineP,failure)
931 : END DO ! iatom_in_halo
932 64 : DEALLOCATE (mat_blk)
933 :
934 : ! clean up
935 64 : DEALLOCATE (atomic_H)
936 64 : DEALLOCATE (atomic_H_blk_row_start)
937 64 : DEALLOCATE (atomic_S)
938 64 : DEALLOCATE (atomic_S_blk_row_start)
939 64 : DEALLOCATE (atomic_filter_mat)
940 :
941 64 : CALL timestop(handle)
942 :
943 384 : END SUBROUTINE fb_fltrmat_add_blkcol_2
944 :
945 : ! **************************************************************************************************
946 : !> \brief generate the list of blocks (atom pairs) to be sent and received
947 : !> in order to construct the filter matrix for each atomic halo.
948 : !> This version is for use with fb_fltrmat_build, where MPI
949 : !> communications are done at each step
950 : !> \param filter_mat : DBCSR formatted filter matrix
951 : !> \param atomic_halo : the halo that contributes to a blk
952 : !> col of the filter matrix
953 : !> \param para_env : cp2k parallel environment
954 : !> \param atom_pairs_send : list of blocks to be sent
955 : !> \param atom_pairs_recv : list of blocks to be received
956 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
957 : ! **************************************************************************************************
958 256 : SUBROUTINE fb_fltrmat_generate_com_pairs(filter_mat, &
959 : atomic_halo, &
960 : para_env, &
961 : atom_pairs_send, &
962 : atom_pairs_recv)
963 : TYPE(dbcsr_type), POINTER :: filter_mat
964 : TYPE(fb_atomic_halo_obj), INTENT(IN) :: atomic_halo
965 : TYPE(mp_para_env_type), POINTER :: para_env
966 : TYPE(fb_com_atom_pairs_obj), INTENT(INOUT) :: atom_pairs_send, atom_pairs_recv
967 :
968 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_generate_com_pairs'
969 :
970 : INTEGER :: dest, handle, iatom_global, &
971 : iatom_in_halo, itask, jatom_global, &
972 : natoms_in_halo, nblkrows_total, &
973 : ntasks_send
974 256 : INTEGER(KIND=int_8), DIMENSION(:, :), POINTER :: tasks_send
975 256 : INTEGER, DIMENSION(:), POINTER :: halo_atoms
976 : TYPE(fb_com_tasks_obj) :: com_tasks_recv, com_tasks_send
977 :
978 256 : CALL timeset(routineN, handle)
979 :
980 256 : NULLIFY (tasks_send)
981 256 : CALL fb_com_tasks_nullify(com_tasks_send)
982 256 : CALL fb_com_tasks_nullify(com_tasks_recv)
983 :
984 : ! initialise atom_pairs_send and atom_pairs_recv
985 256 : IF (fb_com_atom_pairs_has_data(atom_pairs_send)) THEN
986 256 : CALL fb_com_atom_pairs_init(atom_pairs_send)
987 : ELSE
988 0 : CALL fb_com_atom_pairs_create(atom_pairs_send)
989 : END IF
990 256 : IF (fb_com_atom_pairs_has_data(atom_pairs_recv)) THEN
991 256 : CALL fb_com_atom_pairs_init(atom_pairs_recv)
992 : ELSE
993 0 : CALL fb_com_atom_pairs_create(atom_pairs_recv)
994 : END IF
995 :
996 : ! The total number of filter matrix blocks each processor is going
997 : ! to construct equals to the total number of halo atoms in all of
998 : ! the atomic halos local to the processor. The number of send
999 : ! tasks will not exceed this. We do one halo (col) at a time, and
1000 : ! each call of this subroutine will only work on one filter matrix
1001 : ! col corresponding to atomic_halo.
1002 :
1003 : ! The col atom block index for each filter matrix block are the
1004 : ! owner atom of each halo. The row atom block index for each
1005 : ! filter matrix block corresponding to each col are the halo atoms
1006 : ! of the corresponding halos. Filter matrix is non-symmetric: it
1007 : ! is non-square, because the blocks themselves are non-square
1008 :
1009 : CALL fb_atomic_halo_get(atomic_halo=atomic_halo, &
1010 : owner_atom=jatom_global, &
1011 : natoms=natoms_in_halo, &
1012 256 : halo_atoms=halo_atoms)
1013 256 : ntasks_send = natoms_in_halo
1014 :
1015 : ! allocate send tasks
1016 768 : ALLOCATE (tasks_send(TASK_N_RECORDS, ntasks_send))
1017 :
1018 : ! Get the total number of atoms, this can be obtained from the
1019 : ! total number of block rows in the DBCSR filter matrix. We
1020 : ! assumes that before calling this subroutine, the filter_mat has
1021 : ! already been created and initialised: i.e. using
1022 : ! dbcsr_create_new. Even if the matrix is at the moment empty,
1023 : ! the attribute nblkrows_total is already assigned from the dbcsr
1024 : ! distribution data
1025 : CALL dbcsr_get_info(filter_mat, &
1026 256 : nblkrows_total=nblkrows_total)
1027 :
1028 : ! source is always the local processor
1029 : ASSOCIATE (src => para_env%mepos)
1030 : ! construct send tasks
1031 256 : itask = 1
1032 2304 : DO iatom_in_halo = 1, natoms_in_halo
1033 2048 : iatom_global = halo_atoms(iatom_in_halo)
1034 : ! find where the constructed block of filter matrix belongs to
1035 : CALL dbcsr_get_stored_coordinates(filter_mat, &
1036 : iatom_global, &
1037 : jatom_global, &
1038 2048 : processor=dest)
1039 : ! create the send tasks
1040 2048 : tasks_send(TASK_DEST, itask) = dest
1041 2048 : tasks_send(TASK_SRC, itask) = src
1042 : CALL fb_com_tasks_encode_pair(tasks_send(TASK_PAIR, itask), &
1043 : iatom_global, jatom_global, &
1044 2048 : nblkrows_total)
1045 : ! calculation of cost not implemented at the moment
1046 2048 : tasks_send(TASK_COST, itask) = 0
1047 4352 : itask = itask + 1
1048 : END DO ! iatom_in_halo
1049 : END ASSOCIATE
1050 :
1051 256 : CALL fb_com_tasks_create(com_tasks_recv)
1052 256 : CALL fb_com_tasks_create(com_tasks_send)
1053 :
1054 : CALL fb_com_tasks_set(com_tasks=com_tasks_send, &
1055 : task_dim=TASK_N_RECORDS, &
1056 : ntasks=ntasks_send, &
1057 : nencode=nblkrows_total, &
1058 256 : tasks=tasks_send)
1059 :
1060 : ! generate the recv task list (tasks_recv) from the send task list
1061 : CALL fb_com_tasks_transpose_dest_src(com_tasks_recv, "<", com_tasks_send, &
1062 256 : para_env)
1063 :
1064 : ! task lists are now complete, now construct the atom_pairs_send
1065 : ! and atom_pairs_recv from the tasks lists
1066 : CALL fb_com_tasks_build_atom_pairs(com_tasks=com_tasks_send, &
1067 : atom_pairs=atom_pairs_send, &
1068 : natoms_encode=nblkrows_total, &
1069 256 : send_or_recv="send")
1070 : CALL fb_com_tasks_build_atom_pairs(com_tasks=com_tasks_recv, &
1071 : atom_pairs=atom_pairs_recv, &
1072 : natoms_encode=nblkrows_total, &
1073 256 : send_or_recv="recv")
1074 :
1075 : ! cleanup
1076 256 : CALL fb_com_tasks_release(com_tasks_recv)
1077 256 : CALL fb_com_tasks_release(com_tasks_send)
1078 :
1079 256 : CALL timestop(handle)
1080 :
1081 768 : END SUBROUTINE fb_fltrmat_generate_com_pairs
1082 :
1083 : ! **************************************************************************************************
1084 : !> \brief generate the list of blocks (atom pairs) to be sent and received
1085 : !> in order to construct the filter matrix for each atomic halo.
1086 : !> This version is for use with fb_fltrmat_build_2, where MPI
1087 : !> communications are done collectively.
1088 : !> \param filter_mat : DBCSR formatted filter matrix
1089 : !> \param atomic_halos : set of all local atomic halos contributing to the
1090 : !> filter matrix
1091 : !> \param para_env : cp2k parallel environment
1092 : !> \param atom_pairs_send : list of blocks to be sent
1093 : !> \param atom_pairs_recv : list of blocks to be received
1094 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
1095 : ! **************************************************************************************************
1096 16 : SUBROUTINE fb_fltrmat_generate_com_pairs_2(filter_mat, &
1097 : atomic_halos, &
1098 : para_env, &
1099 : atom_pairs_send, &
1100 : atom_pairs_recv)
1101 : TYPE(dbcsr_type), POINTER :: filter_mat
1102 : TYPE(fb_atomic_halo_list_obj), INTENT(IN) :: atomic_halos
1103 : TYPE(mp_para_env_type), POINTER :: para_env
1104 : TYPE(fb_com_atom_pairs_obj), INTENT(INOUT) :: atom_pairs_send, atom_pairs_recv
1105 :
1106 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_generate_com_pairs_2'
1107 :
1108 : INTEGER :: dest, handle, iatom_global, iatom_in_halo, iatom_stored, ihalo, itask, &
1109 : jatom_global, jatom_stored, natoms_in_halo, nblkrows_total, nhalos, ntasks_send
1110 16 : INTEGER(KIND=int_8), DIMENSION(:, :), POINTER :: tasks_send
1111 16 : INTEGER, DIMENSION(:), POINTER :: halo_atoms
1112 : LOGICAL :: transpose
1113 16 : TYPE(fb_atomic_halo_obj), DIMENSION(:), POINTER :: halos
1114 : TYPE(fb_com_tasks_obj) :: com_tasks_recv, com_tasks_send
1115 :
1116 16 : CALL timeset(routineN, handle)
1117 :
1118 16 : NULLIFY (tasks_send)
1119 16 : CALL fb_com_tasks_nullify(com_tasks_send)
1120 16 : CALL fb_com_tasks_nullify(com_tasks_recv)
1121 :
1122 : ! initialise atom_pairs_send and atom_pairs_recv
1123 16 : IF (fb_com_atom_pairs_has_data(atom_pairs_send)) THEN
1124 16 : CALL fb_com_atom_pairs_init(atom_pairs_send)
1125 : ELSE
1126 0 : CALL fb_com_atom_pairs_create(atom_pairs_send)
1127 : END IF
1128 16 : IF (fb_com_atom_pairs_has_data(atom_pairs_recv)) THEN
1129 16 : CALL fb_com_atom_pairs_init(atom_pairs_recv)
1130 : ELSE
1131 0 : CALL fb_com_atom_pairs_create(atom_pairs_recv)
1132 : END IF
1133 :
1134 : ! The col atom block index for each filter matrix block are the
1135 : ! owner atom of each halo. The row atom block index for each
1136 : ! filter matrix block corresponding to each col are the halo atoms
1137 : ! of the corresponding halos. Filter matrix is non-symmetric: it
1138 : ! is non-square, because the blocks themselves are non-square
1139 :
1140 : CALL fb_atomic_halo_list_get(atomic_halos=atomic_halos, &
1141 : nhalos=nhalos, &
1142 16 : halos=halos)
1143 :
1144 : ! estimate the maximum number of blocks (i.e. atom paris) to send
1145 16 : ntasks_send = 0
1146 80 : DO ihalo = 1, nhalos
1147 : CALL fb_atomic_halo_get(atomic_halo=halos(ihalo), &
1148 64 : natoms=natoms_in_halo)
1149 80 : ntasks_send = ntasks_send + natoms_in_halo
1150 : END DO ! ihalo
1151 :
1152 : ! allocate send tasks
1153 48 : ALLOCATE (tasks_send(TASK_N_RECORDS, ntasks_send))
1154 :
1155 : ! Get the total number of atoms. This can be obtained from the
1156 : ! total number of block rows in the DBCSR filter matrix. We
1157 : ! assumes that before calling this subroutine, the filter_mat has
1158 : ! already been created and initialised: i.e. using
1159 : ! dbcsr_create_new. Even if the matrix is at the moment empty,
1160 : ! the attribute nblkrows_total is already assigned from the dbcsr
1161 : ! distribution data
1162 : CALL dbcsr_get_info(filter_mat, &
1163 16 : nblkrows_total=nblkrows_total)
1164 :
1165 : ! source is always the local processor
1166 : ASSOCIATE (src => para_env%mepos)
1167 : ! construct send tasks
1168 16 : itask = 1
1169 80 : DO ihalo = 1, nhalos
1170 : CALL fb_atomic_halo_get(atomic_halo=halos(ihalo), &
1171 : owner_atom=jatom_global, &
1172 : natoms=natoms_in_halo, &
1173 64 : halo_atoms=halo_atoms)
1174 592 : DO iatom_in_halo = 1, natoms_in_halo
1175 512 : iatom_global = halo_atoms(iatom_in_halo)
1176 512 : iatom_stored = iatom_global
1177 512 : jatom_stored = jatom_global
1178 512 : transpose = .FALSE.
1179 : ! find where the constructed block of filter matrix belongs to
1180 : CALL dbcsr_get_stored_coordinates(filter_mat, &
1181 : iatom_stored, &
1182 : jatom_stored, &
1183 512 : processor=dest)
1184 : ! create the send tasks
1185 512 : tasks_send(TASK_DEST, itask) = dest
1186 512 : tasks_send(TASK_SRC, itask) = src
1187 : CALL fb_com_tasks_encode_pair(tasks_send(TASK_PAIR, itask), &
1188 : iatom_global, jatom_global, &
1189 512 : nblkrows_total)
1190 : ! calculation of cost not implemented at the moment
1191 512 : tasks_send(TASK_COST, itask) = 0
1192 1088 : itask = itask + 1
1193 : END DO ! iatom_in_halo
1194 : END DO ! ihalo
1195 : END ASSOCIATE
1196 :
1197 : ! get the actual number of tasks
1198 16 : ntasks_send = itask - 1
1199 :
1200 16 : CALL fb_com_tasks_create(com_tasks_send)
1201 : CALL fb_com_tasks_set(com_tasks=com_tasks_send, &
1202 : task_dim=TASK_N_RECORDS, &
1203 : ntasks=ntasks_send, &
1204 : nencode=nblkrows_total, &
1205 16 : tasks=tasks_send)
1206 :
1207 : ! generate the recv task list (tasks_recv) from the send task list
1208 16 : CALL fb_com_tasks_create(com_tasks_recv)
1209 : CALL fb_com_tasks_transpose_dest_src(com_tasks_recv, "<", com_tasks_send, &
1210 16 : para_env)
1211 :
1212 : ! task lists are now complete, now construct the atom_pairs_send
1213 : ! and atom_pairs_recv from the tasks lists
1214 : CALL fb_com_tasks_build_atom_pairs(com_tasks=com_tasks_send, &
1215 : atom_pairs=atom_pairs_send, &
1216 : natoms_encode=nblkrows_total, &
1217 16 : send_or_recv="send")
1218 : CALL fb_com_tasks_build_atom_pairs(com_tasks=com_tasks_recv, &
1219 : atom_pairs=atom_pairs_recv, &
1220 : natoms_encode=nblkrows_total, &
1221 16 : send_or_recv="recv")
1222 :
1223 : ! cleanup
1224 16 : CALL fb_com_tasks_release(com_tasks_recv)
1225 16 : CALL fb_com_tasks_release(com_tasks_send)
1226 :
1227 16 : CALL timestop(handle)
1228 :
1229 48 : END SUBROUTINE fb_fltrmat_generate_com_pairs_2
1230 :
1231 : ! **************************************************************************************************
1232 : !> \brief Build the atomic filter matrix for each atomic halo
1233 : !> \param atomic_H : atomic KS matrix
1234 : !> \param atomic_S : atomic overlap matrix
1235 : !> \param fermi_level : fermi level used to construct the Fermi-Dirac
1236 : !> filter function
1237 : !> \param filter_temp : temperature used to construct the Fermi-Dirac
1238 : !> filter function
1239 : !> \param atomic_filter_mat : the atomic filter matrix
1240 : !> \param tolerance : anything smaller than tolerance is treated as zero
1241 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
1242 : ! **************************************************************************************************
1243 960 : SUBROUTINE fb_fltrmat_build_atomic_fltrmat(atomic_H, &
1244 320 : atomic_S, &
1245 : fermi_level, &
1246 : filter_temp, &
1247 320 : atomic_filter_mat, &
1248 : tolerance)
1249 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN) :: atomic_H, atomic_S
1250 : REAL(KIND=dp), INTENT(IN) :: fermi_level, filter_temp
1251 : REAL(KIND=dp), DIMENSION(:, :), INTENT(OUT) :: atomic_filter_mat
1252 : REAL(KIND=dp), INTENT(IN) :: tolerance
1253 :
1254 : CHARACTER(LEN=*), PARAMETER :: routineN = 'fb_fltrmat_build_atomic_fltrmat'
1255 :
1256 : CHARACTER(LEN=default_string_length) :: err
1257 : INTEGER :: handle, handle_dgemm, handle_dsygv, ii, &
1258 : info, jj, mat_dim, work_array_size
1259 : LOGICAL :: check_ok
1260 320 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues, filter_function, work
1261 320 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: atomic_S_copy, eigenvectors, &
1262 320 : filtered_eigenvectors
1263 :
1264 320 : CALL timeset(routineN, handle)
1265 :
1266 : ! This subroutine assumes atomic_filter_mat is not zero size, in
1267 : ! other words, it really has to be constructed, instead of just
1268 : ! being a dummy
1269 :
1270 : check_ok = SIZE(atomic_filter_mat, 1) > 0 .AND. &
1271 320 : SIZE(atomic_filter_mat, 2) > 0
1272 0 : CPASSERT(check_ok)
1273 :
1274 : ! initialise
1275 3494720 : atomic_filter_mat = 0.0_dp
1276 320 : mat_dim = SIZE(atomic_H, 1)
1277 :
1278 : ! diagonalise using LAPACK
1279 960 : ALLOCATE (eigenvalues(mat_dim))
1280 : ! get optimal work array size
1281 320 : ALLOCATE (work(1))
1282 : ! dsygv will overwrite part of atomic_H and atomic_S, thus need to copy them
1283 1280 : ALLOCATE (atomic_S_copy(SIZE(atomic_S, 1), SIZE(atomic_S, 2)))
1284 3494720 : atomic_S_copy(:, :) = atomic_S(:, :)
1285 1280 : ALLOCATE (eigenvectors(SIZE(atomic_H, 1), SIZE(atomic_H, 2)))
1286 3494720 : eigenvectors(:, :) = atomic_H(:, :)
1287 :
1288 320 : CALL timeset("fb_atomic_filter_dsygv", handle_dsygv)
1289 :
1290 320 : info = 0
1291 : CALL dsygv(1, 'V', 'U', &
1292 : mat_dim, eigenvectors, mat_dim, &
1293 : atomic_S_copy, mat_dim, eigenvalues, &
1294 320 : work, -1, info)
1295 320 : work_array_size = NINT(work(1))
1296 : ! now allocate work array
1297 320 : DEALLOCATE (work)
1298 960 : ALLOCATE (work(work_array_size))
1299 320 : work = 0.0_dp
1300 : ! do calculation
1301 3494720 : atomic_S_copy(:, :) = atomic_S(:, :)
1302 3494720 : eigenvectors(:, :) = atomic_H(:, :)
1303 320 : info = 0
1304 : CALL dsygv(1, 'V', 'U', &
1305 : mat_dim, eigenvectors, mat_dim, &
1306 : atomic_S_copy, mat_dim, eigenvalues, &
1307 320 : work, work_array_size, info)
1308 : ! check if diagonalisation is successful
1309 320 : IF (info /= 0) THEN
1310 0 : WRITE (err, *) "DSYGV failed with non-zero exit status ", info
1311 0 : CPABORT(TRIM(err))
1312 : END IF
1313 :
1314 320 : CALL timestop(handle_dsygv)
1315 :
1316 320 : DEALLOCATE (work)
1317 320 : DEALLOCATE (atomic_S_copy)
1318 :
1319 : ! first get the filter function
1320 960 : ALLOCATE (filter_function(mat_dim))
1321 320 : filter_function = 0.0_dp
1322 : CALL fb_fltrmat_fermi_dirac_mu(filter_function, &
1323 : eigenvalues, &
1324 : filter_temp, &
1325 320 : fermi_level)
1326 320 : DEALLOCATE (eigenvalues)
1327 :
1328 : ! atomic_H has the eigenvectors, construct the version of it
1329 : ! filtered through the filter function
1330 1280 : ALLOCATE (filtered_eigenvectors(mat_dim, mat_dim))
1331 33600 : DO jj = 1, mat_dim
1332 3494720 : DO ii = 1, mat_dim
1333 : filtered_eigenvectors(ii, jj) = &
1334 3494400 : filter_function(jj)*eigenvectors(ii, jj)
1335 : END DO ! ii
1336 : END DO ! jj
1337 :
1338 320 : DEALLOCATE (filter_function)
1339 :
1340 320 : CALL timeset("fb_atomic_filter_dgemm", handle_dgemm)
1341 :
1342 : ! construct atomic filter matrix
1343 : CALL dgemm("N", &
1344 : "T", &
1345 : mat_dim, &
1346 : mat_dim, &
1347 : mat_dim, &
1348 : 1.0_dp, &
1349 : filtered_eigenvectors, &
1350 : mat_dim, &
1351 : eigenvectors, &
1352 : mat_dim, &
1353 : 0.0_dp, &
1354 : atomic_filter_mat, &
1355 320 : mat_dim)
1356 :
1357 320 : CALL timestop(handle_dgemm)
1358 :
1359 : ! remove small negative terms due to numerical error, the filter
1360 : ! matrix must not be negative definite
1361 33600 : DO jj = 1, SIZE(atomic_filter_mat, 2)
1362 3494720 : DO ii = 1, SIZE(atomic_filter_mat, 1)
1363 3494400 : IF (ABS(atomic_filter_mat(ii, jj)) < tolerance) THEN
1364 52688 : atomic_filter_mat(ii, jj) = 0.0_dp
1365 : END IF
1366 : END DO
1367 : END DO
1368 :
1369 320 : DEALLOCATE (filtered_eigenvectors)
1370 320 : DEALLOCATE (eigenvectors)
1371 :
1372 320 : CALL timestop(handle)
1373 :
1374 960 : END SUBROUTINE fb_fltrmat_build_atomic_fltrmat
1375 :
1376 : ! **************************************************************************************************
1377 : !> \brief get values of Fermi-Dirac distribution based on a given fermi
1378 : !> level at a given set of energy eigenvalues
1379 : !> \param f : the Fermi-Dirac distribution function values
1380 : !> \param eigenvals : set of energy eigenvalues
1381 : !> \param T : temperature
1382 : !> \param mu : the fermi level
1383 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
1384 : ! **************************************************************************************************
1385 320 : SUBROUTINE fb_fltrmat_fermi_dirac_mu(f, eigenvals, T, mu)
1386 : REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: f
1387 : REAL(KIND=dp), DIMENSION(:), INTENT(IN) :: eigenvals
1388 : REAL(KIND=dp), INTENT(IN) :: T, mu
1389 :
1390 : REAL(KIND=dp) :: kTS, ne
1391 :
1392 : ! we want fermi function max at 1, so maxocc = 1 here
1393 :
1394 320 : CALL SmearOcc(f, ne, kTS, eigenvals, mu, T, 1.0_dp, smear_fermi_dirac)
1395 320 : END SUBROUTINE fb_fltrmat_fermi_dirac_mu
1396 :
1397 : ! **************************************************************************************************
1398 : !> \brief get values of Fermi-Dirac distribution based on a given electron
1399 : !> number at a given set of energy eigenvales
1400 : !> \param f : the Fermi-Dirac distribution function values
1401 : !> \param eigenvals : set of energy eigenvalues
1402 : !> \param T : temperature
1403 : !> \param ne : number of electrons
1404 : !> \param maxocc : maximum occupancy per orbital
1405 : !> \author Lianheng Tong (LT) lianheng.tong@kcl.ac.uk
1406 : ! **************************************************************************************************
1407 0 : SUBROUTINE fb_fltrmat_fermi_dirac_ne(f, eigenvals, T, ne, maxocc)
1408 : REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: f
1409 : REAL(KIND=dp), DIMENSION(:), INTENT(IN) :: eigenvals
1410 : REAL(KIND=dp), INTENT(IN) :: T, ne, maxocc
1411 :
1412 : REAL(KIND=dp) :: kTS, mu
1413 :
1414 : ! mu is the calculated fermi level
1415 : ! kTS is the calculated entropic contribution to the energy i.e. -TS
1416 : ! kTS = kT*[f ln f + (1-f) ln (1-f)]
1417 :
1418 0 : CALL SmearFixed(f, mu, kTS, eigenvals, ne, T, maxocc, smear_fermi_dirac)
1419 0 : END SUBROUTINE fb_fltrmat_fermi_dirac_ne
1420 :
1421 : END MODULE qs_fb_filter_matrix_methods
|