(git:281339c)
Loading...
Searching...
No Matches
cp_cfm_dlaf_api.F
Go to the documentation of this file.
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
9
11 USE cp_cfm_types, ONLY: cp_cfm_type
12#if defined(__DLAF)
14 USE dlaf_fortran, ONLY: dlaf_pzheevd, &
15 dlaf_pzhegvd, &
16 dlaf_pzpotrf
17#endif
18 USE kinds, ONLY: dp
19#include "../base/base_uses.f90"
20
21 IMPLICIT NONE
22
23 PRIVATE
24
25 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_cfm_dlaf_api'
26
27 PUBLIC :: cp_cfm_pzpotrf_dlaf
29
30CONTAINS
31
32!***************************************************************************************************
33!> \brief Cholesky factorization using DLA-Future
34!> \param uplo ...
35!> \param n Matrix size
36!> \param a Local matrix
37!> \param ia Row index of first row (has to be 1)
38!> \param ja Col index of first column (has to be 1)
39!> \param desca ScaLAPACK matrix descriptor
40!> \param info 0 if factorization completed normally
41!> \author Rocco Meli
42! **************************************************************************************************
43 SUBROUTINE cp_cfm_pzpotrf_dlaf(uplo, n, a, ia, ja, desca, info)
44 CHARACTER, INTENT(IN) :: uplo
45 INTEGER, INTENT(IN) :: n
46 COMPLEX(KIND=dp), DIMENSION(:, :), TARGET :: a
47 INTEGER, INTENT(IN) :: ia, ja
48 INTEGER, DIMENSION(9) :: desca
49 INTEGER, TARGET :: info
50
51 CHARACTER(len=*), PARAMETER :: routinen = 'cp_cfm_pzpotrf_dlaf'
52
53 INTEGER :: handle
54
55 CALL timeset(routinen, handle)
56#if defined(__DLAF)
57 CALL dlaf_pzpotrf(uplo, n, a, ia, ja, desca, info)
58#else
59 mark_used(uplo)
60 mark_used(n)
61 mark_used(a)
62 mark_used(ia)
63 mark_used(ja)
64 mark_used(desca)
65 mark_used(info)
66 cpabort("CP2K compiled without the DLA-Future library.")
67#endif
68 CALL timestop(handle)
69 END SUBROUTINE cp_cfm_pzpotrf_dlaf
70
71 ! **************************************************************************************************
72!> \brief DLA-Future eigensolver for complex Hermitian matrices
73!> \param matrix ...
74!> \param eigenvectors ...
75!> \param eigenvalues ...
76!> \author Rocco Meli
77! **************************************************************************************************
78 SUBROUTINE cp_cfm_diag_dlaf(matrix, eigenvectors, eigenvalues)
79
80 TYPE(cp_cfm_type), INTENT(IN) :: matrix, eigenvectors
81 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: eigenvalues
82
83 CHARACTER(LEN=*), PARAMETER :: routinen = 'cp_cfm_diag_dlaf'
84
85 INTEGER :: handle, n, nmo
86 REAL(kind=dp), ALLOCATABLE, DIMENSION(:), TARGET :: eig
87
88 CALL timeset(routinen, handle)
89
90 n = matrix%matrix_struct%nrow_global
91 ALLOCATE (eig(n))
92
93 CALL cp_cfm_diag_dlaf_base(matrix, eigenvectors, eig)
94
95 nmo = SIZE(eigenvalues, 1)
96 IF (nmo > n) THEN
97 eigenvalues(1:n) = eig(1:n)
98 ELSE
99 eigenvalues(1:nmo) = eig(1:nmo)
100 END IF
101
102 DEALLOCATE (eig)
103
104 CALL timestop(handle)
105
106 END SUBROUTINE cp_cfm_diag_dlaf
107
108! **************************************************************************************************
109!> \brief DLA-Future generalized eigensolver for complex Hermitian matrices
110!> \param amatrix ...
111!> \param bmatrix ...
112!> \param eigenvectors ...
113!> \param eigenvalues ...
114!> \author Rocco Meli
115! **************************************************************************************************
116 SUBROUTINE cp_cfm_diag_gen_dlaf(amatrix, bmatrix, eigenvectors, eigenvalues)
117
118 TYPE(cp_cfm_type), INTENT(IN) :: amatrix, bmatrix, eigenvectors
119 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: eigenvalues
120
121 CHARACTER(LEN=*), PARAMETER :: routinen = 'cp_cfm_diag_gen_dlaf'
122
123 INTEGER :: handle, n, nmo
124 REAL(kind=dp), ALLOCATABLE, DIMENSION(:), TARGET :: eig
125
126 CALL timeset(routinen, handle)
127
128 n = amatrix%matrix_struct%nrow_global
129 ALLOCATE (eig(n))
130
131 CALL cp_cfm_diag_gen_dlaf_base(amatrix, bmatrix, eigenvectors, eig)
132
133 nmo = SIZE(eigenvalues, 1)
134 IF (nmo > n) THEN
135 eigenvalues(1:n) = eig(1:n)
136 ELSE
137 eigenvalues(1:nmo) = eig(1:nmo)
138 END IF
139
140 DEALLOCATE (eig)
141
142 CALL timestop(handle)
143
144 END SUBROUTINE cp_cfm_diag_gen_dlaf
145
146 !***************************************************************************************************
147!> \brief DLA-Future standard eigensolver for complex Hermitian matrices
148!> \param matrix ...
149!> \param eigenvectors ...
150!> \param eigenvalues ...
151!> \author Rocco Meli
152! **************************************************************************************************
153 SUBROUTINE cp_cfm_diag_dlaf_base(matrix, eigenvectors, eigenvalues)
154 TYPE(cp_cfm_type), INTENT(IN) :: matrix, eigenvectors
155 REAL(kind=dp), DIMENSION(:), INTENT(OUT), TARGET :: eigenvalues
156
157 CHARACTER(len=*), PARAMETER :: dlaf_name = 'pzheevd_dlaf', routinen = 'cp_cfm_diag_dlaf_base'
158 CHARACTER, PARAMETER :: uplo = 'L'
159
160 CHARACTER(LEN=100) :: message
161 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER :: a, z
162 INTEGER :: blacs_context, dlaf_handle, handle, n
163 INTEGER, DIMENSION(9) :: desca, descz
164 INTEGER, TARGET :: info
165
166 CALL timeset(routinen, handle)
167
168#if defined(__DLAF)
169 ! DLAF needs the lower triangular part
170 ! Use eigenvectors matrix as workspace
171 CALL cp_cfm_uplo_to_full(matrix, eigenvectors)
172
173 blacs_context = matrix%matrix_struct%context%get_handle()
174 CALL cp_dlaf_create_grid(blacs_context)
175
176 n = matrix%matrix_struct%nrow_global
177
178 a => matrix%local_data
179 z => eigenvectors%local_data
180
181 desca(:) = matrix%matrix_struct%descriptor(:)
182 descz(:) = eigenvectors%matrix_struct%descriptor(:)
183
184 info = -1
185 CALL timeset(dlaf_name, dlaf_handle)
186 CALL dlaf_pzheevd(uplo, n, a, 1, 1, desca, eigenvalues, z, 1, 1, descz, info)
187 CALL timestop(dlaf_handle)
188
189 IF (info /= 0) THEN
190 WRITE (message, "(A,I0,A)") "ERROR in DLAF_PZHEEVD: Eigensolver failed (INFO = ", info, ")"
191 cpabort(trim(message))
192 END IF
193#else
194 mark_used(a)
195 mark_used(z)
196 mark_used(desca)
197 mark_used(descz)
198 mark_used(matrix)
199 mark_used(eigenvectors)
200 mark_used(eigenvalues)
201 mark_used(uplo)
202 mark_used(n)
203 mark_used(info)
204 mark_used(dlaf_handle)
205 mark_used(dlaf_name)
206 mark_used(message)
207 mark_used(blacs_context)
208 cpabort("CP2K compiled without DLAF library.")
209#endif
210
211 CALL timestop(handle)
212
213 END SUBROUTINE cp_cfm_diag_dlaf_base
214
215!***************************************************************************************************
216!> \brief DLA-Future generalized eigensolver for complex Hermitian matrices
217!> \param amatrix ...
218!> \param bmatrix ...
219!> \param eigenvectors ...
220!> \param eigenvalues ...
221!> \author Rocco Meli
222! **************************************************************************************************
223 SUBROUTINE cp_cfm_diag_gen_dlaf_base(amatrix, bmatrix, eigenvectors, eigenvalues)
224 TYPE(cp_cfm_type), INTENT(IN) :: amatrix, bmatrix, eigenvectors
225 REAL(kind=dp), DIMENSION(:), INTENT(OUT), TARGET :: eigenvalues
226
227 CHARACTER(len=*), PARAMETER :: dlaf_name = 'pzhegvd_dlaf', &
228 routinen = 'cp_cfm_diag_gen_dlaf_base'
229 CHARACTER, PARAMETER :: uplo = 'L'
230
231 CHARACTER(LEN=100) :: message
232 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER :: a, b, z
233 INTEGER :: blacs_context, dlaf_handle, handle, n
234 INTEGER, DIMENSION(9) :: desca, descb, descz
235 INTEGER, TARGET :: info
236
237 CALL timeset(routinen, handle)
238
239#if defined(__DLAF)
240 ! DLAF needs the lower triangular part
241 ! Use eigenvectors matrix as workspace
242 CALL cp_cfm_uplo_to_full(amatrix, eigenvectors)
243 CALL cp_cfm_uplo_to_full(bmatrix, eigenvectors)
244
245 blacs_context = amatrix%matrix_struct%context%get_handle()
246 CALL cp_dlaf_create_grid(blacs_context)
247
248 n = amatrix%matrix_struct%nrow_global
249
250 a => amatrix%local_data
251 b => bmatrix%local_data
252 z => eigenvectors%local_data
253
254 desca(:) = amatrix%matrix_struct%descriptor(:)
255 descb(:) = bmatrix%matrix_struct%descriptor(:)
256 descz(:) = eigenvectors%matrix_struct%descriptor(:)
257
258 info = -1
259 CALL timeset(dlaf_name, dlaf_handle)
260 CALL dlaf_pzhegvd(uplo, n, a, 1, 1, desca, b, 1, 1, descb, eigenvalues, z, 1, 1, descz, info)
261 CALL timestop(dlaf_handle)
262
263 IF (info /= 0) THEN
264 WRITE (message, "(A,I0,A)") "ERROR in DLAF_PZHEGVD: Eigensolver failed (INFO = ", info, ")"
265 cpabort(trim(message))
266 END IF
267#else
268 mark_used(a)
269 mark_used(b)
270 mark_used(z)
271 mark_used(desca)
272 mark_used(descb)
273 mark_used(descz)
274 mark_used(amatrix)
275 mark_used(bmatrix)
276 mark_used(eigenvectors)
277 mark_used(eigenvalues)
278 mark_used(uplo)
279 mark_used(n)
280 mark_used(info)
281 mark_used(dlaf_handle)
282 mark_used(dlaf_name)
283 mark_used(message)
284 mark_used(blacs_context)
285 cpabort("CP2K compiled without DLAF library.")
286#endif
287
288 CALL timestop(handle)
289
290 END SUBROUTINE cp_cfm_diag_gen_dlaf_base
291
292END MODULE cp_cfm_dlaf_api
Basic linear algebra operations for complex full matrices.
subroutine, public cp_cfm_uplo_to_full(matrix, workspace, uplo)
...
subroutine, public cp_cfm_diag_dlaf(matrix, eigenvectors, eigenvalues)
DLA-Future eigensolver for complex Hermitian matrices.
subroutine, public cp_cfm_pzpotrf_dlaf(uplo, n, a, ia, ja, desca, info)
Cholesky factorization using DLA-Future.
subroutine, public cp_cfm_diag_gen_dlaf(amatrix, bmatrix, eigenvectors, eigenvalues)
DLA-Future generalized eigensolver for complex Hermitian matrices.
Represents a complex full matrix distributed on many processors.
subroutine, public cp_dlaf_create_grid(blacs_context)
Create DLA-Future grid from BLACS context.
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Represent a complex full matrix.