Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Interface to (sca)lapack for the Cholesky based procedures
10 : !> \author VW
11 : !> \date 2009-11-09
12 : !> \version 0.8
13 : !>
14 : !> <b>Modification history:</b>
15 : !> - Created 2009-11-09
16 : ! **************************************************************************************************
17 : MODULE cp_dbcsr_diag
18 :
19 : USE cp_blacs_env, ONLY: cp_blacs_env_type
20 : USE cp_cfm_diag, ONLY: cp_cfm_heevd
21 : USE cp_cfm_types, ONLY: cp_cfm_create,&
22 : cp_cfm_release,&
23 : cp_cfm_to_fm,&
24 : cp_cfm_type,&
25 : cp_fm_to_cfm
26 : USE cp_dbcsr_api, ONLY: dbcsr_get_info,&
27 : dbcsr_type
28 : USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
29 : copy_fm_to_dbcsr
30 : USE cp_fm_diag, ONLY: choose_eigv_solver,&
31 : cp_fm_power,&
32 : cp_fm_syevx
33 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
34 : cp_fm_struct_release,&
35 : cp_fm_struct_type
36 : USE cp_fm_types, ONLY: cp_fm_create,&
37 : cp_fm_release,&
38 : cp_fm_type
39 : USE kinds, ONLY: dp
40 : USE message_passing, ONLY: mp_para_env_type
41 : #include "base/base_uses.f90"
42 :
43 : IMPLICIT NONE
44 :
45 : PRIVATE
46 :
47 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_diag'
48 :
49 : ! Public subroutines
50 :
51 : PUBLIC :: cp_dbcsr_syevd, &
52 : cp_dbcsr_syevx, &
53 : cp_dbcsr_heevd, &
54 : cp_dbcsr_power
55 :
56 : CONTAINS
57 :
58 : ! **************************************************************************************************
59 : !> \brief ...
60 : !> \param matrix ...
61 : !> \param eigenvectors ...
62 : !> \param eigenvalues ...
63 : !> \param para_env ...
64 : !> \param blacs_env ...
65 : ! **************************************************************************************************
66 51043 : SUBROUTINE cp_dbcsr_syevd(matrix, eigenvectors, eigenvalues, para_env, blacs_env)
67 :
68 : ! Computes all eigenvalues and vectors of a real symmetric matrix
69 : ! should be quite a bit faster than syevx for that case
70 : ! especially in parallel with tightly clustered evals
71 : ! needs more workspace in the worst case, but much better distributed
72 :
73 : TYPE(dbcsr_type) :: matrix, eigenvectors
74 : REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: eigenvalues
75 : TYPE(mp_para_env_type), POINTER :: para_env
76 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
77 :
78 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_dbcsr_syevd'
79 :
80 : INTEGER :: handle, nfullrows_total
81 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
82 : TYPE(cp_fm_type) :: fm_eigenvectors, fm_matrix
83 :
84 51043 : CALL timeset(routineN, handle)
85 :
86 51043 : NULLIFY (fm_struct)
87 51043 : CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total)
88 :
89 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
90 51043 : ncol_global=nfullrows_total, para_env=para_env)
91 51043 : CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
92 51043 : CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
93 51043 : CALL cp_fm_struct_release(fm_struct)
94 :
95 51043 : CALL copy_dbcsr_to_fm(matrix, fm_matrix)
96 :
97 51043 : CALL choose_eigv_solver(fm_matrix, fm_eigenvectors, eigenvalues)
98 :
99 51043 : CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
100 :
101 51043 : CALL cp_fm_release(fm_matrix)
102 51043 : CALL cp_fm_release(fm_eigenvectors)
103 :
104 51043 : CALL timestop(handle)
105 :
106 51043 : END SUBROUTINE cp_dbcsr_syevd
107 :
108 : ! **************************************************************************************************
109 : !> \brief compute eigenvalues and optionally eigenvectors of a real symmetric matrix using scalapack.
110 : !> If eigenvectors are required this routine will replicate a full matrix on each CPU...
111 : !> if more than a handful of vectors are needed, use cp_dbcsr_syevd instead
112 : !> \param matrix ...
113 : !> \param eigenvectors ...
114 : !> \param eigenvalues ...
115 : !> \param neig ...
116 : !> \param work_syevx ...
117 : !> \param para_env ...
118 : !> \param blacs_env ...
119 : !> \par matrix is supposed to be in upper triangular form, and overwritten by this routine
120 : !> neig is the number of vectors needed (default all)
121 : !> work_syevx evec calculation only, is the fraction of the working buffer allowed (1.0 use full buffer)
122 : !> reducing this saves time, but might cause the routine to fail
123 : ! **************************************************************************************************
124 0 : SUBROUTINE cp_dbcsr_syevx(matrix, eigenvectors, eigenvalues, neig, work_syevx, &
125 : para_env, blacs_env)
126 :
127 : ! Diagonalise the symmetric n by n matrix using the LAPACK library.
128 :
129 : TYPE(dbcsr_type), POINTER :: matrix
130 : TYPE(dbcsr_type), OPTIONAL, POINTER :: eigenvectors
131 : REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: eigenvalues
132 : INTEGER, INTENT(IN), OPTIONAL :: neig
133 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: work_syevx
134 : TYPE(mp_para_env_type), POINTER :: para_env
135 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
136 :
137 : CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_dbcsr_syevx'
138 :
139 : INTEGER :: handle, n, neig_local
140 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
141 : TYPE(cp_fm_type) :: fm_eigenvectors, fm_matrix
142 :
143 0 : CALL timeset(routineN, handle)
144 :
145 : ! by default all
146 0 : CALL dbcsr_get_info(matrix, nfullrows_total=n)
147 0 : neig_local = n
148 0 : IF (PRESENT(neig)) neig_local = neig
149 0 : IF (neig_local == 0) RETURN
150 :
151 0 : NULLIFY (fm_struct)
152 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=n, &
153 0 : ncol_global=n, para_env=para_env)
154 0 : CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
155 :
156 0 : CALL copy_dbcsr_to_fm(matrix, fm_matrix)
157 :
158 0 : IF (PRESENT(eigenvectors)) THEN
159 0 : CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
160 0 : CALL cp_fm_syevx(fm_matrix, fm_eigenvectors, eigenvalues, neig, work_syevx)
161 0 : CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
162 0 : CALL cp_fm_release(fm_eigenvectors)
163 : ELSE
164 0 : CALL cp_fm_syevx(fm_matrix, eigenvalues=eigenvalues, neig=neig, work_syevx=work_syevx)
165 : END IF
166 :
167 0 : CALL cp_fm_struct_release(fm_struct)
168 0 : CALL cp_fm_release(fm_matrix)
169 :
170 0 : CALL timestop(handle)
171 :
172 0 : END SUBROUTINE cp_dbcsr_syevx
173 :
174 : ! **************************************************************************************************
175 : !> \brief ...
176 : !> \param matrix_re ...
177 : !> \param matrix_im ...
178 : !> \param eigenvectors_re ...
179 : !> \param eigenvectors_im ...
180 : !> \param eigenvalues ...
181 : !> \param para_env ...
182 : !> \param blacs_env ...
183 : ! **************************************************************************************************
184 12616 : SUBROUTINE cp_dbcsr_heevd(matrix_re, matrix_im, eigenvectors_re, eigenvectors_im, &
185 3154 : eigenvalues, para_env, blacs_env)
186 :
187 : TYPE(dbcsr_type), OPTIONAL :: matrix_re, matrix_im, eigenvectors_re, &
188 : eigenvectors_im
189 : REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: eigenvalues
190 : TYPE(mp_para_env_type), POINTER :: para_env
191 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
192 :
193 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_dbcsr_heevd'
194 :
195 : INTEGER :: handle, nfullrows_total
196 : TYPE(cp_cfm_type) :: cfm_eigenvectors, cfm_matrix
197 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
198 : TYPE(cp_fm_type) :: fm_eigenvectors_im, fm_eigenvectors_re, &
199 : fm_matrix_im, fm_matrix_re
200 :
201 3154 : CALL timeset(routineN, handle)
202 :
203 : ! Create full matrix structure.
204 3154 : NULLIFY (fm_struct)
205 3154 : IF (PRESENT(matrix_re)) THEN
206 0 : CALL dbcsr_get_info(matrix_re, nfullrows_total=nfullrows_total)
207 3154 : ELSE IF (PRESENT(matrix_im)) THEN
208 3154 : CALL dbcsr_get_info(matrix_im, nfullrows_total=nfullrows_total)
209 : ELSE
210 0 : CPABORT("Neither matrix_re nor matrix_im are present.")
211 : END IF
212 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
213 3154 : ncol_global=nfullrows_total, para_env=para_env)
214 :
215 : ! Create full real matrices.
216 3154 : IF (PRESENT(matrix_re)) THEN
217 0 : CALL cp_fm_create(fm_matrix_re, fm_struct, name="fm_matrix_re")
218 0 : CALL copy_dbcsr_to_fm(matrix_re, fm_matrix_re)
219 : END IF
220 3154 : IF (PRESENT(matrix_im)) THEN
221 3154 : CALL cp_fm_create(fm_matrix_im, fm_struct, name="fm_matrix_im")
222 3154 : CALL copy_dbcsr_to_fm(matrix_im, fm_matrix_im)
223 : END IF
224 :
225 : ! Combine the two real matrices into a complex matrix.
226 3154 : CALL cp_cfm_create(cfm_matrix, fm_struct, name="cfm_matrix")
227 3154 : IF (PRESENT(matrix_re) .AND. PRESENT(matrix_im)) THEN
228 0 : CALL cp_fm_to_cfm(msourcer=fm_matrix_re, msourcei=fm_matrix_im, mtarget=cfm_matrix)
229 3154 : ELSE IF (PRESENT(matrix_re) .AND. .NOT. PRESENT(matrix_im)) THEN
230 0 : CALL cp_fm_to_cfm(msourcer=fm_matrix_re, mtarget=cfm_matrix)
231 3154 : ELSE IF (.NOT. PRESENT(matrix_re) .AND. PRESENT(matrix_im)) THEN
232 3154 : CALL cp_fm_to_cfm(msourcei=fm_matrix_im, mtarget=cfm_matrix)
233 : ELSE
234 0 : CPABORT("Neither matrix_re nor matrix_im are present.")
235 : END IF
236 3154 : IF (PRESENT(matrix_re)) CALL cp_fm_release(fm_matrix_re)
237 3154 : IF (PRESENT(matrix_im)) CALL cp_fm_release(fm_matrix_im)
238 :
239 : ! Diagnonalize the full complex matrix.
240 3154 : CALL cp_cfm_create(cfm_eigenvectors, fm_struct, name="cfm_eigenvectors")
241 3154 : CALL cp_cfm_heevd(cfm_matrix, cfm_eigenvectors, eigenvalues)
242 3154 : CALL cp_cfm_release(cfm_matrix)
243 :
244 : ! Copy the complex eigenvectors back into two real DBCSR matrices.
245 3154 : IF (PRESENT(eigenvectors_re)) THEN
246 3154 : CALL cp_fm_create(fm_eigenvectors_re, fm_struct, name="fm_eigenvectors_re")
247 3154 : CALL cp_cfm_to_fm(msource=cfm_eigenvectors, mtargetr=fm_eigenvectors_re)
248 3154 : CALL copy_fm_to_dbcsr(fm_eigenvectors_re, eigenvectors_re)
249 3154 : CALL cp_fm_release(fm_eigenvectors_re)
250 : END IF
251 3154 : IF (PRESENT(eigenvectors_im)) THEN
252 3154 : CALL cp_fm_create(fm_eigenvectors_im, fm_struct, name="fm_eigenvectors_im")
253 3154 : CALL cp_cfm_to_fm(msource=cfm_eigenvectors, mtargeti=fm_eigenvectors_im)
254 3154 : CALL copy_fm_to_dbcsr(fm_eigenvectors_im, eigenvectors_im)
255 3154 : CALL cp_fm_release(fm_eigenvectors_im)
256 : END IF
257 :
258 : ! Clean up.
259 3154 : CALL cp_cfm_release(cfm_eigenvectors)
260 3154 : CALL cp_fm_struct_release(fm_struct)
261 :
262 3154 : CALL timestop(handle)
263 :
264 3154 : END SUBROUTINE cp_dbcsr_heevd
265 :
266 : ! **************************************************************************************************
267 : !> \brief ...
268 : !> \param matrix ...
269 : !> \param exponent ...
270 : !> \param threshold ...
271 : !> \param n_dependent ...
272 : !> \param para_env ...
273 : !> \param blacs_env ...
274 : !> \param verbose ...
275 : !> \param eigenvectors ...
276 : !> \param eigenvalues ...
277 : ! **************************************************************************************************
278 420 : SUBROUTINE cp_dbcsr_power(matrix, exponent, threshold, n_dependent, para_env, blacs_env, verbose, eigenvectors, eigenvalues)
279 : TYPE(dbcsr_type), INTENT(INOUT) :: matrix
280 : REAL(dp), INTENT(IN) :: exponent, threshold
281 : INTEGER, INTENT(OUT) :: n_dependent
282 : TYPE(mp_para_env_type), POINTER :: para_env
283 : TYPE(cp_blacs_env_type), POINTER :: blacs_env
284 : LOGICAL, INTENT(IN), OPTIONAL :: verbose
285 : TYPE(dbcsr_type), INTENT(INOUT), OPTIONAL :: eigenvectors
286 : REAL(KIND=dp), DIMENSION(2), INTENT(OUT), OPTIONAL :: eigenvalues
287 :
288 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_dbcsr_power'
289 :
290 : INTEGER :: handle, nfullrows_total
291 : REAL(KIND=dp), DIMENSION(2) :: eigenvalues_prv
292 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
293 : TYPE(cp_fm_type) :: fm_eigenvectors, fm_matrix
294 :
295 84 : CALL timeset(routineN, handle)
296 :
297 84 : NULLIFY (fm_struct)
298 84 : CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total)
299 :
300 : CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
301 84 : ncol_global=nfullrows_total, para_env=para_env)
302 84 : CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
303 84 : CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
304 84 : CALL cp_fm_struct_release(fm_struct)
305 :
306 84 : CALL copy_dbcsr_to_fm(matrix, fm_matrix)
307 :
308 84 : CALL cp_fm_power(fm_matrix, fm_eigenvectors, exponent, threshold, n_dependent, verbose, eigenvalues_prv)
309 :
310 84 : CALL copy_fm_to_dbcsr(fm_matrix, matrix)
311 84 : CALL cp_fm_release(fm_matrix)
312 :
313 84 : IF (PRESENT(eigenvalues)) eigenvalues(:) = eigenvalues_prv
314 84 : IF (PRESENT(eigenvectors)) CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
315 :
316 84 : CALL cp_fm_release(fm_eigenvectors)
317 :
318 84 : CALL timestop(handle)
319 :
320 84 : END SUBROUTINE cp_dbcsr_power
321 :
322 : END MODULE cp_dbcsr_diag
|
