(git:374b731)
Loading...
Searching...
No Matches
cp_dbcsr_cholesky.F
Go to the documentation of this file.
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 Interface to (sca)lapack for the Cholesky based procedures
10!> \author VW
11!> \date 2009-09-08
12!> \version 0.8
13!>
14!> <b>Modification history:</b>
15!> - Created 2009-09-08
16! **************************************************************************************************
28 USE cp_fm_types, ONLY: cp_fm_create,&
31 USE dbcsr_api, ONLY: dbcsr_get_info,&
32 dbcsr_type
34#include "base/base_uses.f90"
35
36 IMPLICIT NONE
37
38 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_cholesky'
39
42
43 PRIVATE
44
45CONTAINS
46
47! **************************************************************************************************
48!> \brief used to replace a symmetric positive def. matrix M with its cholesky
49!> decomposition U: M = U^T * U, with U upper triangular
50!> \param matrix the matrix to replace with its cholesky decomposition
51!> \param n the number of row (and columns) of the matrix &
52!> (defaults to the min(size(matrix)))
53!> \param para_env ...
54!> \param blacs_env ...
55!> \par History
56!> 05.2002 created [JVdV]
57!> 12.2002 updated, added n optional parm [fawzi]
58!> \author Joost
59! **************************************************************************************************
60 SUBROUTINE cp_dbcsr_cholesky_decompose(matrix, n, para_env, blacs_env)
61 TYPE(dbcsr_type) :: matrix
62 INTEGER, INTENT(in), OPTIONAL :: n
63 TYPE(mp_para_env_type), POINTER :: para_env
64 TYPE(cp_blacs_env_type), POINTER :: blacs_env
65
66 CHARACTER(len=*), PARAMETER :: routinen = 'cp_dbcsr_cholesky_decompose'
67
68 INTEGER :: handle, my_n, nfullcols_total, &
69 nfullrows_total
70 TYPE(cp_fm_struct_type), POINTER :: fm_struct
71 TYPE(cp_fm_type) :: fm_matrix
72
73 CALL timeset(routinen, handle)
74
75 NULLIFY (fm_struct)
76 CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
77
78 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
79 ncol_global=nfullcols_total, para_env=para_env)
80 CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
81 CALL cp_fm_struct_release(fm_struct)
82
83 CALL copy_dbcsr_to_fm(matrix, fm_matrix)
84
85 my_n = min(fm_matrix%matrix_struct%nrow_global, &
86 fm_matrix%matrix_struct%ncol_global)
87 IF (PRESENT(n)) THEN
88 cpassert(n <= my_n)
89 my_n = n
90 END IF
91
92 CALL cp_fm_potrf(fm_matrix, my_n)
93
94 CALL copy_fm_to_dbcsr(fm_matrix, matrix)
95
96 CALL cp_fm_release(fm_matrix)
97
98 CALL timestop(handle)
99
100 END SUBROUTINE cp_dbcsr_cholesky_decompose
101
102! **************************************************************************************************
103!> \brief used to replace the cholesky decomposition by the inverse
104!> \param matrix the matrix to invert (must be an upper triangular matrix)
105!> \param n size of the matrix to invert (defaults to the min(size(matrix)))
106!> \param para_env ...
107!> \param blacs_env ...
108!> \param upper_to_full ...
109!> \par History
110!> 05.2002 created [JVdV]
111!> \author Joost VandeVondele
112! **************************************************************************************************
113 SUBROUTINE cp_dbcsr_cholesky_invert(matrix, n, para_env, blacs_env, upper_to_full)
114 TYPE(dbcsr_type) :: matrix
115 INTEGER, INTENT(in), OPTIONAL :: n
116 TYPE(mp_para_env_type), POINTER :: para_env
117 TYPE(cp_blacs_env_type), POINTER :: blacs_env
118 LOGICAL, INTENT(IN) :: upper_to_full
119
120 CHARACTER(len=*), PARAMETER :: routinen = 'cp_dbcsr_cholesky_invert'
121
122 INTEGER :: handle, my_n, nfullcols_total, &
123 nfullrows_total
124 TYPE(cp_fm_struct_type), POINTER :: fm_struct
125 TYPE(cp_fm_type) :: fm_matrix, fm_matrix_tmp
126
127 CALL timeset(routinen, handle)
128
129 NULLIFY (fm_struct)
130 CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
131
132 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
133 ncol_global=nfullrows_total, para_env=para_env)
134 CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
135 CALL cp_fm_struct_release(fm_struct)
136
137 CALL copy_dbcsr_to_fm(matrix, fm_matrix)
138
139 my_n = min(fm_matrix%matrix_struct%nrow_global, &
140 fm_matrix%matrix_struct%ncol_global)
141 IF (PRESENT(n)) THEN
142 cpassert(n <= my_n)
143 my_n = n
144 END IF
145
146 CALL cp_fm_potri(fm_matrix, my_n)
147
148 IF (upper_to_full) THEN
149 CALL cp_fm_create(fm_matrix_tmp, fm_matrix%matrix_struct, name="fm_matrix_tmp")
150 CALL cp_fm_upper_to_full(fm_matrix, fm_matrix_tmp)
151 CALL cp_fm_release(fm_matrix_tmp)
152 END IF
153
154 CALL copy_fm_to_dbcsr(fm_matrix, matrix)
155
156 CALL cp_fm_release(fm_matrix)
157
158 CALL timestop(handle)
159
160 END SUBROUTINE cp_dbcsr_cholesky_invert
161
162! **************************************************************************************************
163!> \brief ...
164!> \param matrix ...
165!> \param neig ...
166!> \param matrixb ...
167!> \param matrixout ...
168!> \param op ...
169!> \param pos ...
170!> \param transa ...
171!> \param para_env ...
172!> \param blacs_env ...
173! **************************************************************************************************
174 SUBROUTINE cp_dbcsr_cholesky_restore(matrix, neig, matrixb, matrixout, op, pos, transa, &
175 para_env, blacs_env)
176 TYPE(dbcsr_type) :: matrix
177 INTEGER, INTENT(IN) :: neig
178 TYPE(dbcsr_type) :: matrixb, matrixout
179 CHARACTER(LEN=*), INTENT(IN) :: op
180 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pos, transa
181 TYPE(mp_para_env_type), POINTER :: para_env
182 TYPE(cp_blacs_env_type), POINTER :: blacs_env
183
184 CHARACTER(len=*), PARAMETER :: routinen = 'cp_dbcsr_cholesky_restore'
185
186 CHARACTER :: chol_pos, chol_transa
187 INTEGER :: handle, nfullcols_total, nfullrows_total
188 TYPE(cp_fm_struct_type), POINTER :: fm_struct
189 TYPE(cp_fm_type) :: fm_matrix, fm_matrixb, fm_matrixout
190
191 CALL timeset(routinen, handle)
192
193 NULLIFY (fm_struct)
194
195 CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
196 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
197 ncol_global=nfullcols_total, para_env=para_env)
198 CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
199 CALL cp_fm_struct_release(fm_struct)
200
201 CALL dbcsr_get_info(matrixb, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
202 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
203 ncol_global=nfullcols_total, para_env=para_env)
204 CALL cp_fm_create(fm_matrixb, fm_struct, name="fm_matrixb")
205 CALL cp_fm_struct_release(fm_struct)
206
207 CALL dbcsr_get_info(matrixout, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
208 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
209 ncol_global=nfullcols_total, para_env=para_env)
210 CALL cp_fm_create(fm_matrixout, fm_struct, name="fm_matrixout")
211 CALL cp_fm_struct_release(fm_struct)
212
213 CALL copy_dbcsr_to_fm(matrix, fm_matrix)
214 CALL copy_dbcsr_to_fm(matrixb, fm_matrixb)
215 !CALL copy_dbcsr_to_fm(matrixout, fm_matrixout)
216
217 IF (op /= "SOLVE" .AND. op /= "MULTIPLY") &
218 cpabort("wrong argument op")
219
220 IF (PRESENT(pos)) THEN
221 SELECT CASE (pos)
222 CASE ("LEFT")
223 chol_pos = 'L'
224 CASE ("RIGHT")
225 chol_pos = 'R'
226 CASE DEFAULT
227 cpabort("wrong argument pos")
228 END SELECT
229 ELSE
230 chol_pos = 'L'
231 END IF
232
233 chol_transa = 'N'
234 IF (PRESENT(transa)) chol_transa = transa
235
236 IF ((fm_matrix%use_sp .NEQV. fm_matrixb%use_sp) .OR. (fm_matrix%use_sp .NEQV. fm_matrixout%use_sp)) &
237 cpabort("not the same precision")
238
239 CALL cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, chol_pos, chol_transa)
240
241 CALL copy_fm_to_dbcsr(fm_matrixout, matrixout)
242
243 CALL cp_fm_release(fm_matrix)
244 CALL cp_fm_release(fm_matrixb)
245 CALL cp_fm_release(fm_matrixout)
246
247 CALL timestop(handle)
248
249 END SUBROUTINE cp_dbcsr_cholesky_restore
250
251END MODULE cp_dbcsr_cholesky
252
methods related to the blacs parallel environment
Interface to (sca)lapack for the Cholesky based procedures.
subroutine, public cp_dbcsr_cholesky_decompose(matrix, n, para_env, blacs_env)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
subroutine, public cp_dbcsr_cholesky_restore(matrix, neig, matrixb, matrixout, op, pos, transa, para_env, blacs_env)
...
subroutine, public cp_dbcsr_cholesky_invert(matrix, n, para_env, blacs_env, upper_to_full)
used to replace the cholesky decomposition by the inverse
DBCSR operations in CP2K.
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
basic linear algebra operations for full matrices
subroutine, public cp_fm_potrf(fm_matrix, n)
Cholesky decomposition.
subroutine, public cp_fm_upper_to_full(matrix, work)
given an upper triangular matrix computes the corresponding full matrix
subroutine, public cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, pos, transa)
...
subroutine, public cp_fm_potri(fm_matrix, n)
Invert trianguar matrix.
represent the structure of a full matrix
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Interface to the message passing library MPI.
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
keeps the information about the structure of a full matrix
represent a full matrix
stores all the informations relevant to an mpi environment