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 : ! **************************************************************************************************
9 : !> \brief groups fairly general SCF methods, so that modules other than qs_scf can use them too
10 : !> split off from qs_scf to reduce dependencies
11 : !> \par History
12 : !> - Joost VandeVondele (03.2006)
13 : !> - combine_ks_matrices added (05.04.06,MK)
14 : !> - second ROKS scheme added (15.04.06,MK)
15 : !> - MO occupation management moved (29.08.2008,MK)
16 : !> - correct_mo_eigenvalues was moved from qs_mo_types;
17 : !> new subroutine shift_unocc_mos (03.2016, Sergey Chulkov)
18 : ! **************************************************************************************************
19 : MODULE qs_scf_methods
20 : USE cp_dbcsr_api, ONLY: &
21 : dbcsr_add, dbcsr_desymmetrize, dbcsr_get_block_p, dbcsr_iterator_blocks_left, &
22 : dbcsr_iterator_next_block, dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, &
23 : dbcsr_multiply, dbcsr_p_type, dbcsr_type
24 : USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
25 : cp_dbcsr_sm_fm_multiply
26 : USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale,&
27 : cp_fm_symm,&
28 : cp_fm_triangular_multiply,&
29 : cp_fm_uplo_to_full
30 : USE cp_fm_cholesky, ONLY: cp_fm_cholesky_reduce,&
31 : cp_fm_cholesky_restore
32 : USE cp_fm_diag, ONLY: choose_eigv_solver,&
33 : cp_fm_block_jacobi,&
34 : cp_fm_geeig
35 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
36 : cp_fm_struct_equivalent,&
37 : cp_fm_struct_release,&
38 : cp_fm_struct_type
39 : USE cp_fm_types, ONLY: cp_fm_create,&
40 : cp_fm_get_info,&
41 : cp_fm_release,&
42 : cp_fm_to_fm,&
43 : cp_fm_type
44 : USE input_constants, ONLY: cholesky_inverse,&
45 : cholesky_off,&
46 : cholesky_reduce,&
47 : cholesky_restore
48 : USE kinds, ONLY: dp
49 : USE message_passing, ONLY: mp_para_env_type
50 : USE parallel_gemm_api, ONLY: parallel_gemm
51 : USE qs_density_mixing_types, ONLY: mixing_storage_type
52 : USE qs_mo_types, ONLY: get_mo_set,&
53 : mo_set_type
54 : #include "./base/base_uses.f90"
55 :
56 : IMPLICIT NONE
57 :
58 : PRIVATE
59 :
60 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_scf_methods'
61 : REAL(KIND=dp), PARAMETER :: ratio = 0.25_dp
62 :
63 : PUBLIC :: combine_ks_matrices, &
64 : cp_sm_mix, &
65 : eigensolver, &
66 : eigensolver_generalized, &
67 : eigensolver_dbcsr, &
68 : eigensolver_symm, &
69 : eigensolver_simple, &
70 : scf_env_density_mixing
71 :
72 : INTERFACE combine_ks_matrices
73 : MODULE PROCEDURE combine_ks_matrices_1, &
74 : combine_ks_matrices_2
75 : END INTERFACE combine_ks_matrices
76 :
77 : CONTAINS
78 :
79 : ! **************************************************************************************************
80 : !> \brief perform (if requested) a density mixing
81 : !> \param p_mix_new New density matrices
82 : !> \param mixing_store ...
83 : !> \param rho_ao Density environment
84 : !> \param para_env ...
85 : !> \param iter_delta ...
86 : !> \param iter_count ...
87 : !> \param diis ...
88 : !> \param invert Invert mixing
89 : !> \par History
90 : !> 02.2003 created [fawzi]
91 : !> 08.2014 adapted for kpoints [JGH]
92 : !> \author fawzi
93 : ! **************************************************************************************************
94 158867 : SUBROUTINE scf_env_density_mixing(p_mix_new, mixing_store, rho_ao, para_env, &
95 : iter_delta, iter_count, diis, invert)
96 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: p_mix_new
97 : TYPE(mixing_storage_type), POINTER :: mixing_store
98 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao
99 : TYPE(mp_para_env_type), POINTER :: para_env
100 : REAL(KIND=dp), INTENT(INOUT) :: iter_delta
101 : INTEGER, INTENT(IN) :: iter_count
102 : LOGICAL, INTENT(in), OPTIONAL :: diis, invert
103 :
104 : CHARACTER(len=*), PARAMETER :: routineN = 'scf_env_density_mixing'
105 :
106 : INTEGER :: handle, ic, ispin
107 : LOGICAL :: my_diis, my_invert
108 : REAL(KIND=dp) :: my_p_mix, tmp
109 :
110 158867 : CALL timeset(routineN, handle)
111 :
112 158867 : my_diis = .FALSE.
113 158867 : IF (PRESENT(diis)) my_diis = diis
114 158867 : my_invert = .FALSE.
115 158867 : IF (PRESENT(invert)) my_invert = invert
116 158867 : my_p_mix = mixing_store%alpha
117 158867 : IF (my_diis .OR. iter_count < mixing_store%nskip_mixing) THEN
118 109963 : my_p_mix = 1.0_dp
119 : END IF
120 :
121 158867 : iter_delta = 0.0_dp
122 158867 : CPASSERT(ASSOCIATED(p_mix_new))
123 1449320 : DO ic = 1, SIZE(p_mix_new, 2)
124 2896267 : DO ispin = 1, SIZE(p_mix_new, 1)
125 2737400 : IF (my_invert) THEN
126 159956 : CPASSERT(my_p_mix /= 0.0_dp)
127 159956 : IF (my_p_mix /= 1.0_dp) THEN
128 : CALL dbcsr_add(matrix_a=p_mix_new(ispin, ic)%matrix, &
129 : alpha_scalar=1.0_dp/my_p_mix, &
130 : matrix_b=rho_ao(ispin, ic)%matrix, &
131 25078 : beta_scalar=(my_p_mix - 1.0_dp)/my_p_mix)
132 : END IF
133 : ELSE
134 : CALL cp_sm_mix(m1=p_mix_new(ispin, ic)%matrix, &
135 : m2=rho_ao(ispin, ic)%matrix, &
136 : p_mix=my_p_mix, &
137 : delta=tmp, &
138 1286991 : para_env=para_env)
139 1286991 : iter_delta = MAX(iter_delta, tmp)
140 : END IF
141 : END DO
142 : END DO
143 :
144 158867 : CALL timestop(handle)
145 :
146 158867 : END SUBROUTINE scf_env_density_mixing
147 :
148 : ! **************************************************************************************************
149 : !> \brief Diagonalise the Kohn-Sham matrix to get a new set of MO eigen-
150 : !> vectors and MO eigenvalues. ks will be modified
151 : !> \param matrix_ks_fm ...
152 : !> \param mo_set ...
153 : !> \param ortho ...
154 : !> \param work ...
155 : !> \param cholesky_method ...
156 : !> \param do_level_shift activate the level shifting technique
157 : !> \param level_shift amount of shift applied (in a.u.)
158 : !> \param matrix_u_fm matrix U : S (overlap matrix) = U^T * U
159 : !> \param use_jacobi ...
160 : !> \date 01.05.2001
161 : !> \author Matthias Krack
162 : !> \version 1.0
163 : ! **************************************************************************************************
164 233990 : SUBROUTINE eigensolver(matrix_ks_fm, mo_set, ortho, work, &
165 : cholesky_method, do_level_shift, &
166 : level_shift, matrix_u_fm, use_jacobi)
167 : TYPE(cp_fm_type), INTENT(IN) :: matrix_ks_fm
168 : TYPE(mo_set_type), INTENT(IN) :: mo_set
169 : TYPE(cp_fm_type), INTENT(IN) :: ortho, work
170 : INTEGER, INTENT(inout) :: cholesky_method
171 : LOGICAL, INTENT(in) :: do_level_shift
172 : REAL(KIND=dp), INTENT(in) :: level_shift
173 : TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: matrix_u_fm
174 : LOGICAL, INTENT(in) :: use_jacobi
175 :
176 : CHARACTER(len=*), PARAMETER :: routineN = 'eigensolver'
177 :
178 : INTEGER :: handle, homo, nao, nmo
179 116995 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues
180 : TYPE(cp_fm_type), POINTER :: mo_coeff
181 :
182 116995 : CALL timeset(routineN, handle)
183 :
184 116995 : NULLIFY (mo_coeff)
185 116995 : NULLIFY (mo_eigenvalues)
186 :
187 : ! Diagonalise the Kohn-Sham matrix
188 :
189 : CALL get_mo_set(mo_set=mo_set, &
190 : nao=nao, &
191 : nmo=nmo, &
192 : homo=homo, &
193 : eigenvalues=mo_eigenvalues, &
194 116995 : mo_coeff=mo_coeff)
195 :
196 117043 : SELECT CASE (cholesky_method)
197 : CASE (cholesky_reduce)
198 48 : CALL cp_fm_cholesky_reduce(matrix_ks_fm, ortho)
199 :
200 48 : IF (do_level_shift) THEN
201 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks_fm, mo_coeff=mo_coeff, homo=homo, &
202 28 : level_shift=level_shift, is_triangular=.TRUE., matrix_u_fm=matrix_u_fm)
203 : END IF
204 :
205 48 : CALL choose_eigv_solver(matrix_ks_fm, work, mo_eigenvalues)
206 48 : CALL cp_fm_cholesky_restore(work, nmo, ortho, mo_coeff, "SOLVE")
207 48 : IF (do_level_shift) THEN
208 28 : CALL correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
209 : END IF
210 :
211 : CASE (cholesky_restore)
212 115701 : CALL cp_fm_uplo_to_full(matrix_ks_fm, work)
213 : CALL cp_fm_cholesky_restore(matrix_ks_fm, nao, ortho, work, &
214 115701 : "SOLVE", pos="RIGHT")
215 : CALL cp_fm_cholesky_restore(work, nao, ortho, matrix_ks_fm, &
216 115701 : "SOLVE", pos="LEFT", transa="T")
217 :
218 115701 : IF (do_level_shift) THEN
219 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks_fm, mo_coeff=mo_coeff, homo=homo, &
220 68 : level_shift=level_shift, is_triangular=.TRUE., matrix_u_fm=matrix_u_fm)
221 : END IF
222 :
223 115701 : CALL choose_eigv_solver(matrix_ks_fm, work, mo_eigenvalues)
224 115701 : CALL cp_fm_cholesky_restore(work, nmo, ortho, mo_coeff, "SOLVE")
225 :
226 115701 : IF (do_level_shift) THEN
227 68 : CALL correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
228 : END IF
229 :
230 : CASE (cholesky_inverse)
231 1246 : CALL cp_fm_uplo_to_full(matrix_ks_fm, work)
232 :
233 : CALL cp_fm_triangular_multiply(ortho, matrix_ks_fm, side="R", transpose_tr=.FALSE., &
234 1246 : invert_tr=.FALSE., uplo_tr="U", n_rows=nao, n_cols=nao, alpha=1.0_dp)
235 : CALL cp_fm_triangular_multiply(ortho, matrix_ks_fm, side="L", transpose_tr=.TRUE., &
236 1246 : invert_tr=.FALSE., uplo_tr="U", n_rows=nao, n_cols=nao, alpha=1.0_dp)
237 :
238 1246 : IF (do_level_shift) THEN
239 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks_fm, mo_coeff=mo_coeff, homo=homo, &
240 28 : level_shift=level_shift, is_triangular=.TRUE., matrix_u_fm=matrix_u_fm)
241 : END IF
242 :
243 1246 : CALL choose_eigv_solver(matrix_ks_fm, work, mo_eigenvalues)
244 : CALL cp_fm_triangular_multiply(ortho, work, side="L", transpose_tr=.FALSE., &
245 1246 : invert_tr=.FALSE., uplo_tr="U", n_rows=nao, n_cols=nmo, alpha=1.0_dp)
246 1246 : CALL cp_fm_to_fm(work, mo_coeff, nmo, 1, 1)
247 :
248 118241 : IF (do_level_shift) THEN
249 28 : CALL correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
250 : END IF
251 :
252 : END SELECT
253 :
254 116995 : IF (use_jacobi) THEN
255 0 : CALL cp_fm_to_fm(mo_coeff, ortho)
256 0 : cholesky_method = cholesky_off
257 : END IF
258 :
259 116995 : CALL timestop(handle)
260 :
261 116995 : END SUBROUTINE eigensolver
262 :
263 : ! **************************************************************************************************
264 : !> \brief Solve the generalized eigenvalue problem
265 : !> \param matrix_ks_fm Kohn-Sham matrix in FM format
266 : !> \param matrix_s Overlap matrix (DBCSR)
267 : !> \param mo_set Molecular orbital set
268 : !> \param work Work matrix (used as eigenvector buffer)
269 : ! **************************************************************************************************
270 0 : SUBROUTINE eigensolver_generalized(matrix_ks_fm, matrix_s, mo_set, work)
271 : TYPE(cp_fm_type), INTENT(INOUT) :: matrix_ks_fm
272 : TYPE(dbcsr_type), INTENT(IN) :: matrix_s
273 : TYPE(mo_set_type), INTENT(IN) :: mo_set
274 : TYPE(cp_fm_type), INTENT(INOUT) :: work
275 :
276 : CHARACTER(len=*), PARAMETER :: routineN = 'eigensolver_generalized'
277 :
278 : INTEGER :: handle, nmo
279 0 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues
280 : TYPE(cp_fm_struct_type), POINTER :: fm_struct
281 : TYPE(cp_fm_type) :: s_fm
282 : TYPE(cp_fm_type), POINTER :: mo_coeff
283 :
284 0 : CALL timeset(routineN, handle)
285 :
286 0 : NULLIFY (mo_coeff)
287 0 : NULLIFY (mo_eigenvalues)
288 :
289 0 : CALL get_mo_set(mo_set=mo_set, nmo=nmo, eigenvalues=mo_eigenvalues, mo_coeff=mo_coeff)
290 0 : CALL cp_fm_get_info(matrix_ks_fm, matrix_struct=fm_struct)
291 :
292 : ! Convert S matrix from DBCSR to FM (required for cuSOLVERMp)
293 0 : CALL cp_fm_create(s_fm, fm_struct)
294 0 : CALL copy_dbcsr_to_fm(matrix_s, s_fm)
295 :
296 : ! Solve generalized eigenvalue problem - eigenvectors output to work buffer
297 0 : CALL cp_fm_geeig(matrix_ks_fm, s_fm, mo_coeff, mo_eigenvalues, work)
298 :
299 0 : CALL cp_fm_release(s_fm)
300 :
301 0 : CALL timestop(handle)
302 :
303 0 : END SUBROUTINE eigensolver_generalized
304 :
305 : ! **************************************************************************************************
306 : !> \brief ...
307 : !> \param matrix_ks ...
308 : !> \param matrix_ks_fm ...
309 : !> \param mo_set ...
310 : !> \param ortho_dbcsr ...
311 : !> \param ksbuf1 ...
312 : !> \param ksbuf2 ...
313 : ! **************************************************************************************************
314 8504 : SUBROUTINE eigensolver_dbcsr(matrix_ks, matrix_ks_fm, mo_set, ortho_dbcsr, ksbuf1, ksbuf2)
315 : TYPE(dbcsr_type), INTENT(IN) :: matrix_ks
316 : TYPE(cp_fm_type), INTENT(INOUT) :: matrix_ks_fm
317 : TYPE(mo_set_type), INTENT(IN) :: mo_set
318 : TYPE(dbcsr_type), INTENT(IN) :: ortho_dbcsr
319 : TYPE(dbcsr_type), INTENT(INOUT) :: ksbuf1, ksbuf2
320 :
321 : CHARACTER(len=*), PARAMETER :: routineN = 'eigensolver_dbcsr'
322 :
323 : INTEGER :: handle, nao, nmo
324 2126 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues
325 : TYPE(cp_fm_type) :: all_evecs, nmo_evecs
326 : TYPE(cp_fm_type), POINTER :: mo_coeff
327 :
328 2126 : CALL timeset(routineN, handle)
329 :
330 2126 : NULLIFY (mo_coeff)
331 2126 : NULLIFY (mo_eigenvalues)
332 :
333 : CALL get_mo_set(mo_set=mo_set, &
334 : nao=nao, &
335 : nmo=nmo, &
336 : eigenvalues=mo_eigenvalues, &
337 2126 : mo_coeff=mo_coeff)
338 :
339 : ! Reduce KS matrix
340 2126 : CALL dbcsr_desymmetrize(matrix_ks, ksbuf2)
341 2126 : CALL dbcsr_multiply('N', 'N', 1.0_dp, ksbuf2, ortho_dbcsr, 0.0_dp, ksbuf1)
342 2126 : CALL dbcsr_multiply('T', 'N', 1.0_dp, ortho_dbcsr, ksbuf1, 0.0_dp, ksbuf2)
343 :
344 : ! Solve the eigenvalue problem
345 2126 : CALL copy_dbcsr_to_fm(ksbuf2, matrix_ks_fm)
346 2126 : CALL cp_fm_create(all_evecs, matrix_ks_fm%matrix_struct)
347 2126 : CALL choose_eigv_solver(matrix_ks_fm, all_evecs, mo_eigenvalues)
348 :
349 : ! select first nmo eigenvectors
350 2126 : CALL cp_fm_create(nmo_evecs, mo_coeff%matrix_struct)
351 2126 : CALL cp_fm_to_fm(msource=all_evecs, mtarget=nmo_evecs, ncol=nmo)
352 2126 : CALL cp_fm_release(all_evecs)
353 :
354 : ! Restore the eigenvector of the general eig. problem
355 2126 : CALL cp_dbcsr_sm_fm_multiply(ortho_dbcsr, nmo_evecs, mo_coeff, nmo)
356 :
357 2126 : CALL cp_fm_release(nmo_evecs)
358 2126 : CALL timestop(handle)
359 :
360 2126 : END SUBROUTINE eigensolver_dbcsr
361 :
362 : ! **************************************************************************************************
363 : !> \brief ...
364 : !> \param matrix_ks_fm ...
365 : !> \param mo_set ...
366 : !> \param ortho ...
367 : !> \param work ...
368 : !> \param do_level_shift activate the level shifting technique
369 : !> \param level_shift amount of shift applied (in a.u.)
370 : !> \param matrix_u_fm matrix U : S (overlap matrix) = U^T * U
371 : !> \param use_jacobi ...
372 : !> \param jacobi_threshold ...
373 : !> \param ortho_red ...
374 : !> \param work_red ...
375 : !> \param matrix_ks_fm_red ...
376 : !> \param matrix_u_fm_red ...
377 : ! **************************************************************************************************
378 532 : SUBROUTINE eigensolver_symm(matrix_ks_fm, mo_set, ortho, work, do_level_shift, &
379 : level_shift, matrix_u_fm, use_jacobi, jacobi_threshold, &
380 : ortho_red, work_red, matrix_ks_fm_red, matrix_u_fm_red)
381 : TYPE(cp_fm_type), INTENT(IN) :: matrix_ks_fm
382 : TYPE(mo_set_type), INTENT(IN) :: mo_set
383 : TYPE(cp_fm_type), INTENT(IN) :: ortho, work
384 : LOGICAL, INTENT(IN) :: do_level_shift
385 : REAL(KIND=dp), INTENT(IN) :: level_shift
386 : TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: matrix_u_fm
387 : LOGICAL, INTENT(IN) :: use_jacobi
388 : REAL(KIND=dp), INTENT(IN) :: jacobi_threshold
389 : TYPE(cp_fm_type), INTENT(INOUT), OPTIONAL :: ortho_red, work_red, matrix_ks_fm_red, &
390 : matrix_u_fm_red
391 :
392 : CHARACTER(len=*), PARAMETER :: routineN = 'eigensolver_symm'
393 :
394 : INTEGER :: handle, homo, nao, nao_red, nelectron, &
395 : nmo
396 532 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues
397 532 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues
398 : TYPE(cp_fm_type) :: work_red2
399 : TYPE(cp_fm_type), POINTER :: mo_coeff
400 :
401 532 : CALL timeset(routineN, handle)
402 :
403 532 : NULLIFY (mo_coeff)
404 532 : NULLIFY (mo_eigenvalues)
405 :
406 : ! Diagonalise the Kohn-Sham matrix
407 :
408 : CALL get_mo_set(mo_set=mo_set, &
409 : nao=nao, &
410 : nmo=nmo, &
411 : homo=homo, &
412 : nelectron=nelectron, &
413 : eigenvalues=mo_eigenvalues, &
414 532 : mo_coeff=mo_coeff)
415 :
416 532 : IF (use_jacobi) THEN
417 :
418 0 : CALL cp_fm_symm("L", "U", nao, homo, 1.0_dp, matrix_ks_fm, mo_coeff, 0.0_dp, work)
419 : CALL parallel_gemm("T", "N", homo, nao - homo, nao, 1.0_dp, work, mo_coeff, &
420 0 : 0.0_dp, matrix_ks_fm, b_first_col=homo + 1, c_first_col=homo + 1)
421 :
422 : ! Block Jacobi (pseudo-diagonalization, only one sweep)
423 : CALL cp_fm_block_jacobi(matrix_ks_fm, mo_coeff, mo_eigenvalues, &
424 0 : jacobi_threshold, homo + 1)
425 :
426 : ELSE ! full S^(-1/2) has been computed
427 532 : IF (PRESENT(work_red) .AND. PRESENT(ortho_red) .AND. PRESENT(matrix_ks_fm_red)) THEN
428 532 : CALL cp_fm_get_info(ortho_red, ncol_global=nao_red)
429 532 : CALL cp_fm_symm("L", "U", nao, nao_red, 1.0_dp, matrix_ks_fm, ortho_red, 0.0_dp, work_red)
430 532 : CALL parallel_gemm("T", "N", nao_red, nao_red, nao, 1.0_dp, ortho_red, work_red, 0.0_dp, matrix_ks_fm_red)
431 :
432 532 : IF (do_level_shift) THEN
433 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks_fm_red, mo_coeff=mo_coeff, homo=homo, &
434 86 : level_shift=level_shift, is_triangular=.FALSE., matrix_u_fm=matrix_u_fm_red)
435 : END IF
436 :
437 532 : CALL cp_fm_create(work_red2, matrix_ks_fm_red%matrix_struct)
438 1596 : ALLOCATE (eigenvalues(nao_red))
439 532 : CALL choose_eigv_solver(matrix_ks_fm_red, work_red2, eigenvalues)
440 6244 : mo_eigenvalues(1:MIN(nao_red, nmo)) = eigenvalues(1:MIN(nao_red, nmo))
441 : CALL parallel_gemm("N", "N", nao, nmo, nao_red, 1.0_dp, ortho_red, work_red2, 0.0_dp, &
442 532 : mo_coeff)
443 1596 : CALL cp_fm_release(work_red2)
444 : ELSE
445 0 : CALL cp_fm_symm("L", "U", nao, nao, 1.0_dp, matrix_ks_fm, ortho, 0.0_dp, work)
446 0 : CALL parallel_gemm("T", "N", nao, nao, nao, 1.0_dp, ortho, work, 0.0_dp, matrix_ks_fm)
447 0 : IF (do_level_shift) THEN
448 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks_fm, mo_coeff=mo_coeff, homo=homo, &
449 0 : level_shift=level_shift, is_triangular=.FALSE., matrix_u_fm=matrix_u_fm)
450 : END IF
451 0 : CALL choose_eigv_solver(matrix_ks_fm, work, mo_eigenvalues)
452 : CALL parallel_gemm("N", "N", nao, nmo, nao, 1.0_dp, ortho, work, 0.0_dp, &
453 0 : mo_coeff)
454 : END IF
455 :
456 532 : IF (do_level_shift) THEN
457 86 : CALL correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
458 : END IF
459 :
460 : END IF
461 :
462 532 : CALL timestop(handle)
463 :
464 1064 : END SUBROUTINE eigensolver_symm
465 :
466 : ! **************************************************************************************************
467 :
468 : ! **************************************************************************************************
469 : !> \brief ...
470 : !> \param matrix_ks ...
471 : !> \param mo_set ...
472 : !> \param work ...
473 : !> \param do_level_shift activate the level shifting technique
474 : !> \param level_shift amount of shift applied (in a.u.)
475 : !> \param use_jacobi ...
476 : !> \param jacobi_threshold ...
477 : ! **************************************************************************************************
478 37296 : SUBROUTINE eigensolver_simple(matrix_ks, mo_set, work, do_level_shift, &
479 : level_shift, use_jacobi, jacobi_threshold)
480 :
481 : TYPE(cp_fm_type), INTENT(IN) :: matrix_ks
482 : TYPE(mo_set_type), INTENT(IN) :: mo_set
483 : TYPE(cp_fm_type), INTENT(IN) :: work
484 : LOGICAL, INTENT(IN) :: do_level_shift
485 : REAL(KIND=dp), INTENT(IN) :: level_shift
486 : LOGICAL, INTENT(IN) :: use_jacobi
487 : REAL(KIND=dp), INTENT(IN) :: jacobi_threshold
488 :
489 : CHARACTER(len=*), PARAMETER :: routineN = 'eigensolver_simple'
490 :
491 : INTEGER :: handle, homo, nao, nelectron, nmo
492 18648 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues
493 : TYPE(cp_fm_type), POINTER :: mo_coeff
494 :
495 18648 : CALL timeset(routineN, handle)
496 :
497 18648 : NULLIFY (mo_coeff)
498 18648 : NULLIFY (mo_eigenvalues)
499 :
500 : CALL get_mo_set(mo_set=mo_set, &
501 : nao=nao, &
502 : nmo=nmo, &
503 : homo=homo, &
504 : nelectron=nelectron, &
505 : eigenvalues=mo_eigenvalues, &
506 18648 : mo_coeff=mo_coeff)
507 :
508 18648 : IF (do_level_shift) THEN
509 : ! matrix_u_fm is simply an identity matrix, so we omit it here
510 : CALL shift_unocc_mos(matrix_ks_fm=matrix_ks, mo_coeff=mo_coeff, homo=homo, &
511 0 : level_shift=level_shift, is_triangular=.FALSE.)
512 : END IF
513 :
514 18648 : IF (use_jacobi) THEN
515 18 : CALL cp_fm_symm("L", "U", nao, homo, 1.0_dp, matrix_ks, mo_coeff, 0.0_dp, work)
516 : CALL parallel_gemm("T", "N", homo, nao - homo, nao, 1.0_dp, work, mo_coeff, &
517 18 : 0.0_dp, matrix_ks, b_first_col=homo + 1, c_first_col=homo + 1)
518 : ! Block Jacobi (pseudo-diagonalization, only one sweep)
519 18 : CALL cp_fm_block_jacobi(matrix_ks, mo_coeff, mo_eigenvalues, jacobi_threshold, homo + 1)
520 : ELSE
521 :
522 18630 : CALL choose_eigv_solver(matrix_ks, work, mo_eigenvalues)
523 :
524 18630 : CALL cp_fm_to_fm(work, mo_coeff, nmo, 1, 1)
525 :
526 : END IF
527 :
528 18648 : IF (do_level_shift) THEN
529 0 : CALL correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
530 : END IF
531 :
532 18648 : CALL timestop(handle)
533 :
534 18648 : END SUBROUTINE eigensolver_simple
535 :
536 : ! **************************************************************************************************
537 : !> \brief Perform a mixing of the given matrixes into the first matrix
538 : !> m1 = m2 + p_mix (m1-m2)
539 : !> \param m1 first (new) matrix, is modified
540 : !> \param m2 the second (old) matrix
541 : !> \param p_mix how much m1 is conserved (0: none, 1: all)
542 : !> \param delta maximum norm of m1-m2
543 : !> \param para_env ...
544 : !> \param m3 ...
545 : !> \par History
546 : !> 02.2003 rewamped [fawzi]
547 : !> \author fawzi
548 : !> \note
549 : !> if you what to store the result in m2 swap m1 and m2 an use
550 : !> (1-pmix) as pmix
551 : !> para_env should be removed (embedded in matrix)
552 : ! **************************************************************************************************
553 3456886 : SUBROUTINE cp_sm_mix(m1, m2, p_mix, delta, para_env, m3)
554 :
555 : TYPE(dbcsr_type), POINTER :: m1, m2
556 : REAL(KIND=dp), INTENT(IN) :: p_mix
557 : REAL(KIND=dp), INTENT(OUT) :: delta
558 : TYPE(mp_para_env_type), POINTER :: para_env
559 : TYPE(dbcsr_type), OPTIONAL, POINTER :: m3
560 :
561 : CHARACTER(len=*), PARAMETER :: routineN = 'cp_sm_mix'
562 :
563 : INTEGER :: handle, i, iblock_col, iblock_row, j
564 : LOGICAL :: found
565 1728443 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: p_delta_block, p_new_block, p_old_block
566 : TYPE(dbcsr_iterator_type) :: iter
567 :
568 1728443 : CALL timeset(routineN, handle)
569 1728443 : delta = 0.0_dp
570 :
571 1728443 : CALL dbcsr_iterator_start(iter, m1)
572 24611933 : DO WHILE (dbcsr_iterator_blocks_left(iter))
573 22883490 : CALL dbcsr_iterator_next_block(iter, iblock_row, iblock_col, p_new_block)
574 : CALL dbcsr_get_block_p(matrix=m2, row=iblock_row, col=iblock_col, &
575 22883490 : BLOCK=p_old_block, found=found)
576 22883490 : CPASSERT(ASSOCIATED(p_old_block))
577 24611933 : IF (PRESENT(m3)) THEN
578 : CALL dbcsr_get_block_p(matrix=m3, row=iblock_row, col=iblock_col, &
579 6258410 : BLOCK=p_delta_block, found=found)
580 6258410 : CPASSERT(ASSOCIATED(p_delta_block))
581 :
582 40208263 : DO j = 1, SIZE(p_new_block, 2)
583 305085133 : DO i = 1, SIZE(p_new_block, 1)
584 264876870 : p_delta_block(i, j) = p_new_block(i, j) - p_old_block(i, j)
585 298826723 : delta = MAX(delta, ABS(p_delta_block(i, j)))
586 : END DO
587 : END DO
588 : ELSE
589 97378816 : DO j = 1, SIZE(p_new_block, 2)
590 732938248 : DO i = 1, SIZE(p_new_block, 1)
591 635559432 : p_new_block(i, j) = p_new_block(i, j) - p_old_block(i, j)
592 635559432 : delta = MAX(delta, ABS(p_new_block(i, j)))
593 716313168 : p_new_block(i, j) = p_old_block(i, j) + p_mix*p_new_block(i, j)
594 : END DO
595 : END DO
596 : END IF
597 : END DO
598 1728443 : CALL dbcsr_iterator_stop(iter)
599 :
600 1728443 : CALL para_env%max(delta)
601 :
602 1728443 : CALL timestop(handle)
603 :
604 1728443 : END SUBROUTINE cp_sm_mix
605 :
606 : ! **************************************************************************************************
607 : !> \brief ...
608 : !> \param ksa ...
609 : !> \param ksb ...
610 : !> \param occa ...
611 : !> \param occb ...
612 : !> \param roks_parameter ...
613 : ! **************************************************************************************************
614 1024 : SUBROUTINE combine_ks_matrices_1(ksa, ksb, occa, occb, roks_parameter)
615 :
616 : ! Combine the alpha and beta Kohn-Sham matrices during a restricted open
617 : ! Kohn-Sham (ROKS) calculation
618 : ! On input ksa and ksb contain the alpha and beta Kohn-Sham matrices,
619 : ! respectively. occa and occb contain the corresponding MO occupation
620 : ! numbers. On output the combined ROKS operator matrix is returned in ksa.
621 :
622 : ! Literature: - C. C. J. Roothaan, Rev. Mod. Phys. 32, 179 (1960)
623 : ! - M. F. Guest and V. R. Saunders, Mol. Phys. 28(3), 819 (1974)
624 :
625 : TYPE(cp_fm_type), INTENT(IN) :: ksa, ksb
626 : REAL(KIND=dp), DIMENSION(:), INTENT(IN) :: occa, occb
627 : REAL(KIND=dp), DIMENSION(0:2, 0:2, 1:2), &
628 : INTENT(IN) :: roks_parameter
629 :
630 : CHARACTER(LEN=*), PARAMETER :: routineN = 'combine_ks_matrices_1'
631 :
632 : INTEGER :: handle, i, icol_global, icol_local, &
633 : irow_global, irow_local, j, &
634 : ncol_local, nrow_local
635 1024 : INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
636 : LOGICAL :: compatible_matrices
637 : REAL(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
638 1024 : POINTER :: fa, fb
639 : TYPE(cp_fm_struct_type), POINTER :: ksa_struct, ksb_struct
640 :
641 : ! -------------------------------------------------------------------------
642 :
643 1024 : CALL timeset(routineN, handle)
644 :
645 : CALL cp_fm_get_info(matrix=ksa, &
646 : matrix_struct=ksa_struct, &
647 : nrow_local=nrow_local, &
648 : ncol_local=ncol_local, &
649 : row_indices=row_indices, &
650 : col_indices=col_indices, &
651 1024 : local_data=fa)
652 :
653 : CALL cp_fm_get_info(matrix=ksb, &
654 : matrix_struct=ksb_struct, &
655 1024 : local_data=fb)
656 :
657 1024 : compatible_matrices = cp_fm_struct_equivalent(ksa_struct, ksb_struct)
658 1024 : CPASSERT(compatible_matrices)
659 :
660 16398 : IF (SUM(occb) == 0.0_dp) fb = 0.0_dp
661 :
662 16398 : DO icol_local = 1, ncol_local
663 15374 : icol_global = col_indices(icol_local)
664 15374 : j = INT(occa(icol_global)) + INT(occb(icol_global))
665 181549 : DO irow_local = 1, nrow_local
666 165151 : irow_global = row_indices(irow_local)
667 165151 : i = INT(occa(irow_global)) + INT(occb(irow_global))
668 : fa(irow_local, icol_local) = &
669 : roks_parameter(i, j, 1)*fa(irow_local, icol_local) + &
670 180525 : roks_parameter(i, j, 2)*fb(irow_local, icol_local)
671 : END DO
672 : END DO
673 :
674 1024 : CALL timestop(handle)
675 :
676 1024 : END SUBROUTINE combine_ks_matrices_1
677 :
678 : ! **************************************************************************************************
679 : !> \brief ...
680 : !> \param ksa ...
681 : !> \param ksb ...
682 : !> \param occa ...
683 : !> \param occb ...
684 : !> \param f ...
685 : !> \param nalpha ...
686 : !> \param nbeta ...
687 : ! **************************************************************************************************
688 0 : SUBROUTINE combine_ks_matrices_2(ksa, ksb, occa, occb, f, nalpha, nbeta)
689 :
690 : ! Combine the alpha and beta Kohn-Sham matrices during a restricted open
691 : ! Kohn-Sham (ROKS) calculation
692 : ! On input ksa and ksb contain the alpha and beta Kohn-Sham matrices,
693 : ! respectively. occa and occb contain the corresponding MO occupation
694 : ! numbers. On output the combined ROKS operator matrix is returned in ksa.
695 :
696 : ! Literature: - C. C. J. Roothaan, Rev. Mod. Phys. 32, 179 (1960)
697 : ! - M. Filatov and S. Shaik, Chem. Phys. Lett. 288, 689 (1998)
698 :
699 : TYPE(cp_fm_type), INTENT(IN) :: ksa, ksb
700 : REAL(KIND=dp), DIMENSION(:), INTENT(IN) :: occa, occb
701 : REAL(KIND=dp), INTENT(IN) :: f
702 : INTEGER, INTENT(IN) :: nalpha, nbeta
703 :
704 : CHARACTER(LEN=*), PARAMETER :: routineN = 'combine_ks_matrices_2'
705 :
706 : INTEGER :: handle, icol_global, icol_local, &
707 : irow_global, irow_local, ncol_local, &
708 : nrow_local
709 0 : INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
710 : LOGICAL :: compatible_matrices
711 : REAL(KIND=dp) :: beta, t1, t2, ta, tb
712 : REAL(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
713 0 : POINTER :: fa, fb
714 : TYPE(cp_fm_struct_type), POINTER :: ksa_struct, ksb_struct
715 :
716 : ! -------------------------------------------------------------------------
717 :
718 0 : CALL timeset(routineN, handle)
719 :
720 : CALL cp_fm_get_info(matrix=ksa, &
721 : matrix_struct=ksa_struct, &
722 : nrow_local=nrow_local, &
723 : ncol_local=ncol_local, &
724 : row_indices=row_indices, &
725 : col_indices=col_indices, &
726 0 : local_data=fa)
727 :
728 : CALL cp_fm_get_info(matrix=ksb, &
729 : matrix_struct=ksb_struct, &
730 0 : local_data=fb)
731 :
732 0 : compatible_matrices = cp_fm_struct_equivalent(ksa_struct, ksb_struct)
733 0 : CPASSERT(compatible_matrices)
734 :
735 0 : beta = 1.0_dp/(1.0_dp - f)
736 :
737 0 : DO icol_local = 1, ncol_local
738 :
739 0 : icol_global = col_indices(icol_local)
740 :
741 0 : DO irow_local = 1, nrow_local
742 :
743 0 : irow_global = row_indices(irow_local)
744 :
745 0 : t1 = 0.5_dp*(fa(irow_local, icol_local) + fb(irow_local, icol_local))
746 :
747 0 : IF ((0 < irow_global) .AND. (irow_global <= nbeta)) THEN
748 0 : IF ((0 < icol_global) .AND. (icol_global <= nbeta)) THEN
749 : ! closed-closed
750 0 : fa(irow_local, icol_local) = t1
751 0 : ELSE IF ((nbeta < icol_global) .AND. (icol_global <= nalpha)) THEN
752 : ! closed-open
753 0 : ta = 0.5_dp*(f - REAL(occa(icol_global), KIND=dp))/f
754 0 : tb = 0.5_dp*(f - REAL(occb(icol_global), KIND=dp))/f
755 0 : t2 = ta*fa(irow_local, icol_local) + tb*fb(irow_local, icol_local)
756 0 : fa(irow_local, icol_local) = t1 + (beta - 1.0_dp)*t2
757 : ELSE
758 : ! closed-virtual
759 0 : fa(irow_local, icol_local) = t1
760 : END IF
761 0 : ELSE IF ((nbeta < irow_global) .AND. (irow_global <= nalpha)) THEN
762 : IF ((0 < irow_global) .AND. (irow_global <= nbeta)) THEN
763 : ! open-closed
764 : ta = 0.5_dp*(f - REAL(occa(irow_global), KIND=dp))/f
765 : tb = 0.5_dp*(f - REAL(occb(irow_global), KIND=dp))/f
766 : t2 = ta*fa(irow_local, icol_local) + tb*fb(irow_local, icol_local)
767 : fa(irow_local, icol_local) = t1 + (beta - 1.0_dp)*t2
768 0 : ELSE IF ((nbeta < icol_global) .AND. (icol_global <= nalpha)) THEN
769 : ! open-open
770 0 : ta = 0.5_dp*(f - REAL(occa(icol_global), KIND=dp))/f
771 0 : tb = 0.5_dp*(f - REAL(occb(icol_global), KIND=dp))/f
772 0 : t2 = ta*fa(irow_local, icol_local) + tb*fb(irow_local, icol_local)
773 0 : IF (irow_global == icol_global) THEN
774 0 : fa(irow_local, icol_local) = t1 - t2
775 : ELSE
776 0 : fa(irow_local, icol_local) = t1 - 0.5_dp*t2
777 : END IF
778 : ELSE
779 : ! open-virtual
780 0 : ta = 0.5_dp*(f - REAL(occa(irow_global), KIND=dp))/f
781 0 : tb = 0.5_dp*(f - REAL(occb(irow_global), KIND=dp))/f
782 0 : t2 = ta*fa(irow_local, icol_local) + tb*fb(irow_local, icol_local)
783 0 : fa(irow_local, icol_local) = t1 - t2
784 : END IF
785 : ELSE
786 0 : IF ((0 < irow_global) .AND. (irow_global < nbeta)) THEN
787 : ! virtual-closed
788 0 : fa(irow_local, icol_local) = t1
789 0 : ELSE IF ((nbeta < icol_global) .AND. (icol_global <= nalpha)) THEN
790 : ! virtual-open
791 0 : ta = 0.5_dp*(f - REAL(occa(icol_global), KIND=dp))/f
792 0 : tb = 0.5_dp*(f - REAL(occb(icol_global), KIND=dp))/f
793 0 : t2 = ta*fa(irow_local, icol_local) + tb*fb(irow_local, icol_local)
794 0 : fa(irow_local, icol_local) = t1 - t2
795 : ELSE
796 : ! virtual-virtual
797 0 : fa(irow_local, icol_local) = t1
798 : END IF
799 : END IF
800 :
801 : END DO
802 : END DO
803 :
804 0 : CALL timestop(handle)
805 :
806 0 : END SUBROUTINE combine_ks_matrices_2
807 :
808 : ! **************************************************************************************************
809 : !> \brief Correct MO eigenvalues after MO level shifting.
810 : !> \param mo_eigenvalues vector of eigenvalues
811 : !> \param homo index of the highest occupied molecular orbital
812 : !> \param nmo number of molecular orbitals
813 : !> \param level_shift amount of applied level shifting (in a.u.)
814 : !> \date 19.04.2002
815 : !> \par History
816 : !> - correct_mo_eigenvalues added (18.04.02,MK)
817 : !> - moved from module qs_mo_types, revised interface (03.2016, Sergey Chulkov)
818 : !> \author MK
819 : !> \version 1.0
820 : ! **************************************************************************************************
821 210 : PURE SUBROUTINE correct_mo_eigenvalues(mo_eigenvalues, homo, nmo, level_shift)
822 :
823 : REAL(kind=dp), DIMENSION(:), INTENT(inout) :: mo_eigenvalues
824 : INTEGER, INTENT(in) :: homo, nmo
825 : REAL(kind=dp), INTENT(in) :: level_shift
826 :
827 : INTEGER :: imo
828 :
829 5530 : DO imo = homo + 1, nmo
830 5530 : mo_eigenvalues(imo) = mo_eigenvalues(imo) - level_shift
831 : END DO
832 :
833 210 : END SUBROUTINE correct_mo_eigenvalues
834 :
835 : ! **************************************************************************************************
836 : !> \brief Adjust the Kohn-Sham matrix by shifting the orbital energies of all
837 : !> unoccupied molecular orbitals
838 : !> \param matrix_ks_fm transformed Kohn-Sham matrix = U^{-1,T} * KS * U^{-1}
839 : !> \param mo_coeff matrix of molecular orbitals (C)
840 : !> \param homo number of occupied molecular orbitals
841 : !> \param level_shift amount of shift applying (in a.u.)
842 : !> \param is_triangular indicates that matrix_u_fm contains an upper triangular matrix
843 : !> \param matrix_u_fm matrix U: S (overlap matrix) = U^T * U;
844 : !> assume an identity matrix if omitted
845 : !> \par History
846 : !> 03.2016 created [Sergey Chulkov]
847 : ! **************************************************************************************************
848 210 : SUBROUTINE shift_unocc_mos(matrix_ks_fm, mo_coeff, homo, &
849 : level_shift, is_triangular, matrix_u_fm)
850 :
851 : TYPE(cp_fm_type), INTENT(IN) :: matrix_ks_fm, mo_coeff
852 : INTEGER, INTENT(in) :: homo
853 : REAL(kind=dp), INTENT(in) :: level_shift
854 : LOGICAL, INTENT(in) :: is_triangular
855 : TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: matrix_u_fm
856 :
857 : CHARACTER(len=*), PARAMETER :: routineN = 'shift_unocc_mos'
858 :
859 : INTEGER :: handle, nao, nao_red, nmo
860 210 : REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: weights
861 : TYPE(cp_fm_struct_type), POINTER :: ao_mo_fmstruct
862 : TYPE(cp_fm_type) :: u_mo, u_mo_scaled
863 :
864 210 : CALL timeset(routineN, handle)
865 :
866 210 : IF (PRESENT(matrix_u_fm)) THEN
867 210 : CALL cp_fm_get_info(mo_coeff, ncol_global=nmo)
868 210 : CALL cp_fm_get_info(matrix_u_fm, nrow_global=nao_red, ncol_global=nao)
869 : ELSE
870 0 : CALL cp_fm_get_info(mo_coeff, ncol_global=nmo, nrow_global=nao)
871 0 : nao_red = nao
872 : END IF
873 :
874 210 : NULLIFY (ao_mo_fmstruct)
875 : CALL cp_fm_struct_create(ao_mo_fmstruct, nrow_global=nao_red, ncol_global=nmo, &
876 210 : para_env=mo_coeff%matrix_struct%para_env, context=mo_coeff%matrix_struct%context)
877 :
878 210 : CALL cp_fm_create(u_mo, ao_mo_fmstruct)
879 210 : CALL cp_fm_create(u_mo_scaled, ao_mo_fmstruct)
880 :
881 210 : CALL cp_fm_struct_release(ao_mo_fmstruct)
882 :
883 : ! U * C
884 210 : IF (PRESENT(matrix_u_fm)) THEN
885 210 : IF (is_triangular) THEN
886 124 : CALL cp_fm_to_fm(mo_coeff, u_mo)
887 : CALL cp_fm_triangular_multiply(matrix_u_fm, u_mo, side="L", transpose_tr=.FALSE., &
888 124 : invert_tr=.FALSE., uplo_tr="U", n_rows=nao, n_cols=nmo, alpha=1.0_dp)
889 : ELSE
890 86 : CALL parallel_gemm("N", "N", nao_red, nmo, nao, 1.0_dp, matrix_u_fm, mo_coeff, 0.0_dp, u_mo)
891 : END IF
892 : ELSE
893 : ! assume U is an identity matrix
894 0 : CALL cp_fm_to_fm(mo_coeff, u_mo)
895 : END IF
896 :
897 210 : CALL cp_fm_to_fm(u_mo, u_mo_scaled)
898 :
899 : ! up-shift all unoccupied molecular orbitals by the amount of 'level_shift'
900 : ! weight = diag(DELTA) = (0, ... 0, level_shift, ..., level_shift)
901 : ! MO index : 1 .. homo homo+1 ... nmo
902 630 : ALLOCATE (weights(nmo))
903 1120 : weights(1:homo) = 0.0_dp
904 5820 : weights(homo + 1:nmo) = level_shift
905 : ! DELTA * U * C
906 : ! DELTA is a diagonal matrix, so simply scale all the columns of (U * C) by weights(:)
907 210 : CALL cp_fm_column_scale(u_mo_scaled, weights)
908 210 : DEALLOCATE (weights)
909 :
910 : ! NewKS = U^{-1,T} * KS * U^{-1} + (U * C) * DELTA * (U * C)^T
911 210 : CALL parallel_gemm("N", "T", nao_red, nao_red, nmo, 1.0_dp, u_mo, u_mo_scaled, 1.0_dp, matrix_ks_fm)
912 :
913 210 : CALL cp_fm_release(u_mo_scaled)
914 210 : CALL cp_fm_release(u_mo)
915 :
916 210 : CALL timestop(handle)
917 :
918 420 : END SUBROUTINE shift_unocc_mos
919 :
920 : END MODULE qs_scf_methods
|