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 Provides interfaces to LAPACK eigenvalue/SVD routines
10 : !> \par History
11 : !> JGH (26-5-2001): delay D/S C/Z problem to the lapack library call
12 : !> \author APSI
13 : ! **************************************************************************************************
14 : MODULE eigenvalueproblems
15 :
16 : USE kinds, ONLY: dp
17 : #include "../base/base_uses.f90"
18 :
19 : IMPLICIT NONE
20 :
21 : PRIVATE
22 :
23 : PUBLIC :: diagonalise
24 :
25 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'eigenvalueproblems'
26 :
27 : INTERFACE diagonalise
28 : MODULE PROCEDURE diagonalise_ssyev
29 : MODULE PROCEDURE diagonalise_chpev
30 : END INTERFACE
31 :
32 : INTERFACE singular_values
33 : MODULE PROCEDURE cp2k_sgesvd
34 : MODULE PROCEDURE cp2k_cgesvd
35 : END INTERFACE
36 :
37 : CONTAINS
38 :
39 : ! **************************************************************************************************
40 : !> \brief ...
41 : !> \param matrix ...
42 : !> \param mysize ...
43 : !> \param uplo ...
44 : !> \param eigenvalues ...
45 : !> \param eigenvectors ...
46 : ! **************************************************************************************************
47 988 : SUBROUTINE diagonalise_ssyev(matrix, mysize, uplo, eigenvalues, &
48 988 : eigenvectors)
49 :
50 : REAL(KIND=dp), INTENT(IN) :: matrix(:, :)
51 : INTEGER, INTENT(IN) :: mysize
52 : CHARACTER(LEN=1), INTENT(IN) :: uplo
53 : REAL(KIND=dp), INTENT(OUT) :: eigenvalues(:), eigenvectors(:, :)
54 :
55 : CHARACTER, PARAMETER :: jobz = "V"
56 :
57 : INTEGER :: info, lda, lwork
58 1976 : REAL(KIND=dp) :: work(3*mysize - 1)
59 :
60 988 : lda = SIZE(matrix, 1)
61 988 : lwork = 3*mysize - 1
62 :
63 12844 : eigenvectors = matrix
64 :
65 : CALL dsyev(jobz, uplo, mysize, eigenvectors, lda, eigenvalues, &
66 988 : work, lwork, info)
67 988 : IF (info /= 0) THEN
68 0 : CPABORT("Error in diagonalisation")
69 : END IF
70 :
71 988 : END SUBROUTINE diagonalise_ssyev
72 :
73 : ! **************************************************************************************************
74 : !> \brief ...
75 : !> \param matrix ...
76 : !> \param mysize ...
77 : !> \param uplo ...
78 : !> \param eigenvalues ...
79 : !> \param eigenvectors ...
80 : ! **************************************************************************************************
81 0 : SUBROUTINE diagonalise_chpev(matrix, mysize, uplo, eigenvalues, &
82 0 : eigenvectors)
83 :
84 : COMPLEX(KIND=dp), INTENT(INOUT) :: matrix(:)
85 : INTEGER, INTENT(IN) :: mysize
86 : CHARACTER(LEN=1), INTENT(IN) :: uplo
87 : REAL(KIND=dp), INTENT(OUT) :: eigenvalues(:)
88 : COMPLEX(KIND=dp), INTENT(OUT) :: eigenvectors(:, :)
89 :
90 : CHARACTER, PARAMETER :: jobz = "V"
91 :
92 : INTEGER :: info
93 0 : COMPLEX(KIND=dp) :: work(2*mysize - 1)
94 0 : REAL(KIND=dp) :: rwork(3*mysize - 2)
95 :
96 : CALL zhpev(jobz, uplo, mysize, matrix, eigenvalues, &
97 0 : eigenvectors, mysize, work, rwork, info)
98 0 : IF (info /= 0) THEN
99 0 : CPABORT("Error in diagonalisation")
100 : END IF
101 :
102 0 : END SUBROUTINE diagonalise_chpev
103 :
104 : ! **************************************************************************************************
105 : !> \brief ...
106 : !> \param matrix ...
107 : !> \param svalues ...
108 : !> \param mrow ...
109 : !> \param ncol ...
110 : !> \param uvec ...
111 : !> \param vtvec ...
112 : ! **************************************************************************************************
113 0 : SUBROUTINE cp2k_sgesvd(matrix, svalues, mrow, ncol, uvec, vtvec)
114 :
115 : REAL(KIND=dp), INTENT(IN) :: matrix(:, :)
116 : REAL(KIND=dp), INTENT(OUT) :: svalues(:)
117 : INTEGER, INTENT(IN) :: mrow, ncol
118 : REAL(KIND=dp), INTENT(OUT) :: uvec(:, :), vtvec(:, :)
119 :
120 : CHARACTER, PARAMETER :: jobu = "A", jobvt = "A"
121 :
122 : INTEGER :: info, lda, ldu, ldvt, lwork
123 0 : REAL(KIND=dp) :: work(25*(mrow + ncol))
124 :
125 0 : lwork = 25*(mrow + ncol)
126 0 : lda = SIZE(matrix, 1)
127 0 : ldu = SIZE(uvec, 1)
128 0 : ldvt = SIZE(vtvec, 1)
129 :
130 : CALL dgesvd(jobu, jobvt, mrow, ncol, matrix, lda, svalues, &
131 0 : uvec, ldu, vtvec, ldvt, work, lwork, info)
132 0 : IF (info /= 0) THEN
133 0 : CPABORT("Error in singular value decomposition.")
134 : END IF
135 :
136 0 : END SUBROUTINE cp2k_sgesvd
137 :
138 : ! **************************************************************************************************
139 : !> \brief ...
140 : !> \param matrix ...
141 : !> \param svalues ...
142 : !> \param mrow ...
143 : !> \param ncol ...
144 : !> \param uvec ...
145 : !> \param vtvec ...
146 : ! **************************************************************************************************
147 0 : SUBROUTINE cp2k_cgesvd(matrix, svalues, mrow, ncol, uvec, vtvec)
148 :
149 : COMPLEX(KIND=dp), INTENT(IN) :: matrix(:, :)
150 : REAL(KIND=dp), INTENT(OUT) :: svalues(:)
151 : INTEGER, INTENT(IN) :: mrow, ncol
152 : COMPLEX(KIND=dp), INTENT(OUT) :: uvec(:, :), vtvec(:, :)
153 :
154 : CHARACTER, PARAMETER :: jobu = "A", jobvt = "A"
155 :
156 : INTEGER :: info, lda, ldu, ldvt, lwork
157 0 : COMPLEX(KIND=dp) :: work(25*(mrow + ncol))
158 0 : REAL(KIND=dp) :: rwork(25*(mrow + ncol))
159 :
160 0 : lwork = 25*(mrow + ncol)
161 0 : lda = SIZE(matrix, 1)
162 0 : ldu = SIZE(uvec, 1)
163 0 : ldvt = SIZE(vtvec, 1)
164 :
165 : CALL zgesvd(jobu, jobvt, mrow, ncol, matrix, lda, svalues, &
166 0 : uvec, ldu, vtvec, ldvt, work, lwork, rwork, info)
167 0 : IF (info /= 0) THEN
168 0 : CPABORT("Error in singular value decomposition.")
169 : END IF
170 :
171 0 : END SUBROUTINE cp2k_cgesvd
172 :
173 : END MODULE eigenvalueproblems
174 :
|
