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 Some utilities for the construction of
10 : !> the localization environment
11 : !> \author MI (05-2005)
12 : ! **************************************************************************************************
13 : MODULE qs_loc_utils
14 :
15 : USE cell_types, ONLY: cell_type
16 : USE cp_array_utils, ONLY: cp_1d_r_p_type
17 : USE cp_control_types, ONLY: dft_control_type
18 : USE cp_dbcsr_api, ONLY: dbcsr_copy,&
19 : dbcsr_p_type,&
20 : dbcsr_type
21 : USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply
22 : USE cp_files, ONLY: close_file,&
23 : open_file
24 : USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale
25 : USE cp_fm_diag, ONLY: choose_eigv_solver
26 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
27 : cp_fm_struct_release,&
28 : cp_fm_struct_type
29 : USE cp_fm_types, ONLY: &
30 : cp_fm_create, cp_fm_get_info, cp_fm_get_submatrix, cp_fm_release, cp_fm_set_all, &
31 : cp_fm_set_submatrix, cp_fm_to_fm, cp_fm_type, cp_fm_write_unformatted
32 : USE cp_log_handling, ONLY: cp_get_default_logger,&
33 : cp_logger_get_default_io_unit,&
34 : cp_logger_type,&
35 : cp_to_string
36 : USE cp_output_handling, ONLY: cp_p_file,&
37 : cp_print_key_finished_output,&
38 : cp_print_key_generate_filename,&
39 : cp_print_key_should_output,&
40 : cp_print_key_unit_nr
41 : USE distribution_1d_types, ONLY: distribution_1d_type
42 : USE input_constants, ONLY: &
43 : do_loc_crazy, do_loc_direct, do_loc_gapo, do_loc_jacobi, do_loc_l1_norm_sd, do_loc_none, &
44 : do_loc_scdm, energy_loc_range, op_loc_berry, op_loc_boys, op_loc_pipek, state_loc_all, &
45 : state_loc_list, state_loc_mixed, state_loc_none, state_loc_range
46 : USE input_section_types, ONLY: section_vals_get_subs_vals,&
47 : section_vals_type,&
48 : section_vals_val_get
49 : USE kinds, ONLY: default_path_length,&
50 : default_string_length,&
51 : dp
52 : USE mathconstants, ONLY: twopi
53 : USE memory_utilities, ONLY: reallocate
54 : USE message_passing, ONLY: mp_para_env_type
55 : USE parallel_gemm_api, ONLY: parallel_gemm
56 : USE particle_types, ONLY: particle_type
57 : USE qs_environment_types, ONLY: get_qs_env,&
58 : qs_environment_type
59 : USE qs_loc_types, ONLY: get_qs_loc_env,&
60 : localized_wfn_control_create,&
61 : localized_wfn_control_release,&
62 : localized_wfn_control_type,&
63 : qs_loc_env_type,&
64 : set_qs_loc_env
65 : USE qs_localization_methods, ONLY: initialize_weights
66 : USE qs_mo_methods, ONLY: make_mo_eig
67 : USE qs_mo_types, ONLY: get_mo_set,&
68 : mo_set_type
69 : USE qs_operators_ao, ONLY: build_exp_ikr_matrix
70 : USE qs_scf_types, ONLY: ot_method_nr
71 : USE scf_control_types, ONLY: scf_control_type
72 : #include "./base/base_uses.f90"
73 :
74 : IMPLICIT NONE
75 :
76 : PRIVATE
77 :
78 : ! *** Global parameters ***
79 :
80 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_loc_utils'
81 :
82 : ! *** Public ***
83 : PUBLIC :: qs_loc_env_init, loc_write_restart, &
84 : retain_history, qs_loc_init, compute_berry_operator, &
85 : set_loc_centers, set_loc_wfn_lists, qs_loc_control_init
86 :
87 : CONTAINS
88 :
89 : ! **************************************************************************************************
90 : !> \brief copy old mos to new ones, allocating as necessary
91 : !> \param mo_loc_history ...
92 : !> \param mo_loc ...
93 : ! **************************************************************************************************
94 10 : SUBROUTINE retain_history(mo_loc_history, mo_loc)
95 :
96 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: mo_loc_history
97 : TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: mo_loc
98 :
99 : CHARACTER(len=*), PARAMETER :: routineN = 'retain_history'
100 :
101 : INTEGER :: handle, i, ncol_hist, ncol_loc
102 :
103 10 : CALL timeset(routineN, handle)
104 :
105 10 : IF (.NOT. ASSOCIATED(mo_loc_history)) THEN
106 8 : ALLOCATE (mo_loc_history(SIZE(mo_loc)))
107 4 : DO i = 1, SIZE(mo_loc_history)
108 4 : CALL cp_fm_create(mo_loc_history(i), mo_loc(i)%matrix_struct)
109 : END DO
110 : END IF
111 :
112 20 : DO i = 1, SIZE(mo_loc_history)
113 10 : CALL cp_fm_get_info(mo_loc_history(i), ncol_global=ncol_hist)
114 10 : CALL cp_fm_get_info(mo_loc(i), ncol_global=ncol_loc)
115 10 : CPASSERT(ncol_hist == ncol_loc)
116 30 : CALL cp_fm_to_fm(mo_loc(i), mo_loc_history(i))
117 : END DO
118 :
119 10 : CALL timestop(handle)
120 :
121 10 : END SUBROUTINE retain_history
122 :
123 : ! **************************************************************************************************
124 : !> \brief rotate the mo_new, so that the orbitals are as similar
125 : !> as possible to ones in mo_ref.
126 : !> \param mo_new ...
127 : !> \param mo_ref ...
128 : !> \param matrix_S ...
129 : ! **************************************************************************************************
130 8 : SUBROUTINE rotate_state_to_ref(mo_new, mo_ref, matrix_S)
131 :
132 : TYPE(cp_fm_type), INTENT(IN) :: mo_new, mo_ref
133 : TYPE(dbcsr_type), POINTER :: matrix_S
134 :
135 : CHARACTER(len=*), PARAMETER :: routineN = 'rotate_state_to_ref'
136 :
137 : INTEGER :: handle, ncol, ncol_ref, nrow
138 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues
139 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
140 : TYPE(cp_fm_type) :: o1, o2, o3, o4, smo
141 :
142 8 : CALL timeset(routineN, handle)
143 :
144 8 : CALL cp_fm_get_info(mo_new, nrow_global=nrow, ncol_global=ncol)
145 8 : CALL cp_fm_get_info(mo_ref, ncol_global=ncol_ref)
146 8 : CPASSERT(ncol == ncol_ref)
147 :
148 8 : NULLIFY (fm_struct_tmp)
149 8 : CALL cp_fm_create(smo, mo_ref%matrix_struct)
150 :
151 : CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, &
152 : ncol_global=ncol, para_env=mo_new%matrix_struct%para_env, &
153 8 : context=mo_new%matrix_struct%context)
154 8 : CALL cp_fm_create(o1, fm_struct_tmp)
155 8 : CALL cp_fm_create(o2, fm_struct_tmp)
156 8 : CALL cp_fm_create(o3, fm_struct_tmp)
157 8 : CALL cp_fm_create(o4, fm_struct_tmp)
158 8 : CALL cp_fm_struct_release(fm_struct_tmp)
159 :
160 : ! o1 = (mo_new)^T matrix_S mo_ref
161 8 : CALL cp_dbcsr_sm_fm_multiply(matrix_S, mo_ref, smo, ncol)
162 8 : CALL parallel_gemm('T', 'N', ncol, ncol, nrow, 1.0_dp, mo_new, smo, 0.0_dp, o1)
163 :
164 : ! o2 = (o1^T o1)
165 8 : CALL parallel_gemm('T', 'N', ncol, ncol, ncol, 1.0_dp, o1, o1, 0.0_dp, o2)
166 :
167 : ! o2 = (o1^T o1)^-1/2
168 24 : ALLOCATE (eigenvalues(ncol))
169 8 : CALL choose_eigv_solver(o2, o3, eigenvalues)
170 8 : CALL cp_fm_to_fm(o3, o4)
171 72 : eigenvalues(:) = 1.0_dp/SQRT(eigenvalues(:))
172 8 : CALL cp_fm_column_scale(o4, eigenvalues)
173 8 : CALL parallel_gemm('N', 'T', ncol, ncol, ncol, 1.0_dp, o3, o4, 0.0_dp, o2)
174 :
175 : ! o3 = o1 (o1^T o1)^-1/2
176 8 : CALL parallel_gemm('N', 'N', ncol, ncol, ncol, 1.0_dp, o1, o2, 0.0_dp, o3)
177 :
178 : ! mo_new o1 (o1^T o1)^-1/2
179 8 : CALL parallel_gemm('N', 'N', nrow, ncol, ncol, 1.0_dp, mo_new, o3, 0.0_dp, smo)
180 8 : CALL cp_fm_to_fm(smo, mo_new)
181 :
182 : ! XXXXXXX testing
183 : ! CALL parallel_gemm('N','T',ncol,ncol,ncol,1.0_dp,o3,o3,0.0_dp,o1)
184 : ! WRITE(*,*) o1%local_data
185 : ! CALL parallel_gemm('T','N',ncol,ncol,ncol,1.0_dp,o3,o3,0.0_dp,o1)
186 : ! WRITE(*,*) o1%local_data
187 :
188 8 : CALL cp_fm_release(o1)
189 8 : CALL cp_fm_release(o2)
190 8 : CALL cp_fm_release(o3)
191 8 : CALL cp_fm_release(o4)
192 8 : CALL cp_fm_release(smo)
193 :
194 8 : CALL timestop(handle)
195 :
196 32 : END SUBROUTINE rotate_state_to_ref
197 :
198 : ! **************************************************************************************************
199 : !> \brief allocates the data, and initializes the operators
200 : !> \param qs_loc_env new environment for the localization calculations
201 : !> \param localized_wfn_control variables and directives for the localization
202 : !> \param qs_env the qs_env in which the qs_env lives
203 : !> \param myspin ...
204 : !> \param do_localize ...
205 : !> \param loc_coeff ...
206 : !> \param mo_loc_history ...
207 : !> \par History
208 : !> 04.2005 created [MI]
209 : !> \author MI
210 : !> \note
211 : !> similar to the old one, but not quite
212 : ! **************************************************************************************************
213 968 : SUBROUTINE qs_loc_env_init(qs_loc_env, localized_wfn_control, qs_env, myspin, do_localize, &
214 484 : loc_coeff, mo_loc_history)
215 :
216 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
217 : TYPE(localized_wfn_control_type), POINTER :: localized_wfn_control
218 : TYPE(qs_environment_type), POINTER :: qs_env
219 : INTEGER, INTENT(IN), OPTIONAL :: myspin
220 : LOGICAL, INTENT(IN), OPTIONAL :: do_localize
221 : TYPE(cp_fm_type), DIMENSION(:), INTENT(IN), &
222 : OPTIONAL :: loc_coeff
223 : TYPE(cp_fm_type), DIMENSION(:), OPTIONAL, POINTER :: mo_loc_history
224 :
225 : CHARACTER(len=*), PARAMETER :: routineN = 'qs_loc_env_init'
226 :
227 : INTEGER :: dim_op, handle, i, iatom, imo, imoloc, &
228 : ispin, j, l_spin, lb, nao, naosub, &
229 : natoms, nmo, nmosub, nspins, s_spin, ub
230 : LOGICAL :: loc_coeff_spin_resolved
231 : REAL(KIND=dp) :: my_occ, occ_imo
232 484 : REAL(KIND=dp), DIMENSION(:), POINTER :: occupations
233 484 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: vecbuffer
234 : TYPE(cell_type), POINTER :: cell
235 : TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
236 484 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: moloc_coeff
237 : TYPE(cp_fm_type), POINTER :: mat_ptr, mo_coeff
238 484 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
239 : TYPE(distribution_1d_type), POINTER :: local_molecules
240 484 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
241 : TYPE(mp_para_env_type), POINTER :: para_env
242 484 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
243 :
244 484 : CALL timeset(routineN, handle)
245 :
246 484 : NULLIFY (mos, matrix_s, moloc_coeff, particle_set, para_env, cell, &
247 484 : local_molecules, occupations, mat_ptr)
248 484 : IF (PRESENT(do_localize)) qs_loc_env%do_localize = do_localize
249 484 : IF (qs_loc_env%do_localize) THEN
250 : CALL get_qs_env(qs_env=qs_env, matrix_s=matrix_s, cell=cell, &
251 : local_molecules=local_molecules, particle_set=particle_set, &
252 484 : para_env=para_env, mos=mos)
253 484 : nspins = SIZE(mos, 1)
254 484 : loc_coeff_spin_resolved = .FALSE.
255 484 : IF (PRESENT(loc_coeff)) THEN
256 308 : loc_coeff_spin_resolved = nspins*2 == SIZE(loc_coeff)
257 : END IF
258 484 : s_spin = 1
259 484 : l_spin = nspins
260 484 : IF (PRESENT(myspin)) THEN
261 162 : s_spin = myspin
262 162 : l_spin = myspin
263 : END IF
264 484 : IF (loc_coeff_spin_resolved) THEN
265 166 : ALLOCATE (moloc_coeff(s_spin:s_spin + 2*(l_spin - s_spin) + 1))
266 : ELSE
267 1948 : ALLOCATE (moloc_coeff(s_spin:l_spin))
268 : END IF
269 1102 : DO ispin = s_spin, l_spin
270 618 : NULLIFY (tmp_fm_struct, mo_coeff)
271 618 : CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff, nao=nao, nmo=nmo)
272 618 : nmosub = localized_wfn_control%nloc_states(ispin)
273 : CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nao, &
274 618 : ncol_global=nmosub, para_env=para_env, context=mo_coeff%matrix_struct%context)
275 618 : IF (loc_coeff_spin_resolved) THEN
276 44 : CALL cp_fm_create(moloc_coeff(2*ispin - 1), tmp_fm_struct)
277 44 : CALL cp_fm_create(moloc_coeff(2*ispin), tmp_fm_struct)
278 : ELSE
279 574 : CALL cp_fm_create(moloc_coeff(ispin), tmp_fm_struct)
280 : END IF
281 :
282 : CALL cp_fm_get_info(moloc_coeff(ispin), nrow_global=naosub, &
283 618 : ncol_global=nmosub)
284 618 : CPASSERT(nao == naosub)
285 618 : IF ((localized_wfn_control%do_homo) .OR. &
286 : (localized_wfn_control%set_of_states == state_loc_mixed)) THEN
287 606 : CPASSERT(nmo >= nmosub)
288 : ELSE
289 12 : CPASSERT(nao - nmo >= nmosub)
290 : END IF
291 618 : CALL cp_fm_set_all(moloc_coeff(ispin), 0.0_dp)
292 2338 : CALL cp_fm_struct_release(tmp_fm_struct)
293 : END DO ! ispin
294 : ! Copy the submatrix
295 :
296 484 : IF (PRESENT(loc_coeff)) ALLOCATE (mat_ptr)
297 :
298 1102 : DO ispin = s_spin, l_spin
299 : CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff, &
300 618 : occupation_numbers=occupations, nao=nao, nmo=nmo)
301 618 : lb = localized_wfn_control%lu_bound_states(1, ispin)
302 618 : ub = localized_wfn_control%lu_bound_states(2, ispin)
303 :
304 618 : IF (PRESENT(loc_coeff)) THEN
305 428 : mat_ptr = loc_coeff(ispin)
306 : ELSE
307 190 : mat_ptr => mo_coeff
308 : END IF
309 618 : IF ((localized_wfn_control%set_of_states == state_loc_list) .OR. &
310 : (localized_wfn_control%set_of_states == state_loc_mixed)) THEN
311 444 : ALLOCATE (vecbuffer(1, nao))
312 148 : IF (localized_wfn_control%do_homo) THEN
313 134 : my_occ = occupations(localized_wfn_control%loc_states(1, ispin))
314 : END IF
315 148 : nmosub = SIZE(localized_wfn_control%loc_states, 1)
316 148 : CPASSERT(nmosub > 0)
317 148 : imoloc = 0
318 934 : DO i = lb, ub
319 : ! Get the index in the subset
320 786 : imoloc = imoloc + 1
321 : ! Get the index in the full set
322 786 : imo = localized_wfn_control%loc_states(i, ispin)
323 786 : IF (localized_wfn_control%do_homo) THEN
324 652 : occ_imo = occupations(imo)
325 652 : IF (ABS(occ_imo - my_occ) > localized_wfn_control%eps_occ) THEN
326 0 : IF (localized_wfn_control%localization_method /= do_loc_none) THEN
327 : CALL cp_abort(__LOCATION__, &
328 : "States with different occupations "// &
329 0 : "cannot be rotated together")
330 : END IF
331 : END IF
332 : END IF
333 : ! Take the imo vector from the full set and copy in the imoloc vector of the subset
334 : CALL cp_fm_get_submatrix(mat_ptr, vecbuffer, 1, imo, &
335 786 : nao, 1, transpose=.TRUE.)
336 : CALL cp_fm_set_submatrix(moloc_coeff(ispin), vecbuffer, 1, imoloc, &
337 934 : nao, 1, transpose=.TRUE.)
338 : END DO
339 148 : DEALLOCATE (vecbuffer)
340 : ELSE
341 470 : my_occ = occupations(lb)
342 470 : occ_imo = occupations(ub)
343 470 : IF (ABS(occ_imo - my_occ) > localized_wfn_control%eps_occ) THEN
344 0 : IF (localized_wfn_control%localization_method /= do_loc_none) THEN
345 : CALL cp_abort(__LOCATION__, &
346 : "States with different occupations "// &
347 0 : "cannot be rotated together")
348 : END IF
349 : END IF
350 470 : nmosub = localized_wfn_control%nloc_states(ispin)
351 :
352 470 : IF (loc_coeff_spin_resolved) THEN
353 44 : CALL cp_fm_to_fm(loc_coeff(2*ispin - 1), moloc_coeff(2*ispin - 1))
354 44 : CALL cp_fm_to_fm(loc_coeff(2*ispin), moloc_coeff(2*ispin))
355 : ELSE
356 426 : CALL cp_fm_to_fm(mat_ptr, moloc_coeff(ispin), nmosub, lb, 1)
357 : END IF
358 : END IF
359 :
360 : ! we have the mo's to be localized now, see if we can rotate them according to the history
361 : ! only do that if we have a history of course. The history is filled
362 1720 : IF (PRESENT(mo_loc_history)) THEN
363 104 : IF (localized_wfn_control%use_history .AND. ASSOCIATED(mo_loc_history)) THEN
364 : CALL rotate_state_to_ref(moloc_coeff(ispin), &
365 8 : mo_loc_history(ispin), matrix_s(1)%matrix)
366 : END IF
367 : END IF
368 :
369 : END DO
370 :
371 484 : IF (PRESENT(loc_coeff)) DEALLOCATE (mat_ptr)
372 :
373 : CALL set_qs_loc_env(qs_loc_env=qs_loc_env, cell=cell, local_molecules=local_molecules, &
374 : moloc_coeff=moloc_coeff, particle_set=particle_set, para_env=para_env, &
375 484 : localized_wfn_control=localized_wfn_control)
376 :
377 : ! Prepare the operators
378 484 : NULLIFY (tmp_fm_struct, mo_coeff)
379 1452 : nmosub = MAXVAL(localized_wfn_control%nloc_states)
380 484 : CALL get_mo_set(mos(1), mo_coeff=mo_coeff)
381 : CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nmosub, &
382 484 : ncol_global=nmosub, para_env=para_env, context=mo_coeff%matrix_struct%context)
383 :
384 484 : IF (localized_wfn_control%operator_type == op_loc_berry) THEN
385 478 : IF (qs_loc_env%cell%orthorhombic) THEN
386 466 : dim_op = 3
387 : ELSE
388 12 : dim_op = 6
389 : END IF
390 478 : CALL set_qs_loc_env(qs_loc_env=qs_loc_env, dim_op=dim_op)
391 5844 : ALLOCATE (qs_loc_env%op_sm_set(2, dim_op))
392 1948 : DO i = 1, dim_op
393 4888 : DO j = 1, SIZE(qs_loc_env%op_sm_set, 1)
394 2940 : NULLIFY (qs_loc_env%op_sm_set(j, i)%matrix)
395 2940 : ALLOCATE (qs_loc_env%op_sm_set(j, i)%matrix)
396 : CALL dbcsr_copy(qs_loc_env%op_sm_set(j, i)%matrix, matrix_s(1)%matrix, &
397 4410 : name="qs_loc_env%op_sm_"//TRIM(ADJUSTL(cp_to_string(j)))//"-"//TRIM(ADJUSTL(cp_to_string(i))))
398 : END DO
399 : END DO
400 :
401 6 : ELSE IF (localized_wfn_control%operator_type == op_loc_pipek) THEN
402 6 : natoms = SIZE(qs_loc_env%particle_set, 1)
403 96 : ALLOCATE (qs_loc_env%op_fm_set(natoms, 1))
404 6 : CALL set_qs_loc_env(qs_loc_env=qs_loc_env, dim_op=natoms)
405 12 : DO ispin = 1, SIZE(qs_loc_env%op_fm_set, 2)
406 6 : CALL get_mo_set(mos(ispin), nmo=nmo)
407 84 : DO iatom = 1, natoms
408 72 : CALL cp_fm_create(qs_loc_env%op_fm_set(iatom, ispin), tmp_fm_struct)
409 :
410 72 : CALL cp_fm_get_info(qs_loc_env%op_fm_set(iatom, ispin), nrow_global=nmosub)
411 72 : CPASSERT(nmo >= nmosub)
412 150 : CALL cp_fm_set_all(qs_loc_env%op_fm_set(iatom, ispin), 0.0_dp)
413 : END DO ! iatom
414 : END DO ! ispin
415 : ELSE
416 0 : CPABORT("Type of operator not implemented")
417 : END IF
418 484 : CALL cp_fm_struct_release(tmp_fm_struct)
419 :
420 484 : IF (localized_wfn_control%operator_type == op_loc_berry) THEN
421 :
422 478 : CALL initialize_weights(qs_loc_env%cell, qs_loc_env%weights)
423 :
424 478 : CALL get_berry_operator(qs_loc_env, qs_env)
425 :
426 : ELSE IF (localized_wfn_control%operator_type == op_loc_pipek) THEN
427 :
428 : !! here we don't have to do anything
429 : !! CALL get_pipek_mezey_operator ( qs_loc_env, qs_env )
430 :
431 : END IF
432 :
433 484 : qs_loc_env%molecular_states = .FALSE.
434 484 : qs_loc_env%wannier_states = .FALSE.
435 : END IF
436 484 : CALL timestop(handle)
437 :
438 484 : END SUBROUTINE qs_loc_env_init
439 :
440 : ! **************************************************************************************************
441 : !> \brief A wrapper to compute the Berry operator for periodic systems
442 : !> \param qs_loc_env new environment for the localization calculations
443 : !> \param qs_env the qs_env in which the qs_env lives
444 : !> \par History
445 : !> 04.2005 created [MI]
446 : !> 04.2018 modified [RZK, ZL]
447 : !> \author MI
448 : ! **************************************************************************************************
449 478 : SUBROUTINE get_berry_operator(qs_loc_env, qs_env)
450 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
451 : TYPE(qs_environment_type), POINTER :: qs_env
452 :
453 : CHARACTER(len=*), PARAMETER :: routineN = 'get_berry_operator'
454 :
455 : INTEGER :: dim_op, handle
456 : TYPE(cell_type), POINTER :: cell
457 478 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: op_sm_set
458 :
459 478 : CALL timeset(routineN, handle)
460 :
461 478 : NULLIFY (cell, op_sm_set)
462 : CALL get_qs_loc_env(qs_loc_env=qs_loc_env, op_sm_set=op_sm_set, &
463 478 : cell=cell, dim_op=dim_op)
464 478 : CALL compute_berry_operator(qs_env, cell, op_sm_set, dim_op)
465 :
466 478 : CALL timestop(handle)
467 478 : END SUBROUTINE get_berry_operator
468 :
469 : ! **************************************************************************************************
470 : !> \brief Computes the reciprocal-space operators used by Berry localization.
471 : !> The operators contain the contracted AO matrix elements of
472 : !> cos(k*r) and sin(k*r) for reciprocal vectors selected from the cell.
473 : !> This routine preserves the historical localization interface and
474 : !> delegates the AO assembly to qs_operators_ao::build_exp_ikr_matrix.
475 : !> \param qs_env the qs_env in which the qs_env lives
476 : !> \param cell cell used to generate the reciprocal vectors and fold positions
477 : !> \param op_sm_set cosine and sine operator storage
478 : !> \param dim_op number of reciprocal vectors, either 3 or 6
479 : !> \par History
480 : !> 04.2005 created [MI]
481 : !> 04.2018 wrapped old code [RZK, ZL]
482 : !> \author MI
483 : !> \note The supplied cell is also used for coordinate folding and its
484 : !> periodicity is temporarily set to three dimensions during assembly.
485 : ! **************************************************************************************************
486 504 : SUBROUTINE compute_berry_operator(qs_env, cell, op_sm_set, dim_op)
487 : TYPE(qs_environment_type), POINTER :: qs_env
488 : TYPE(cell_type), POINTER :: cell
489 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: op_sm_set
490 : INTEGER :: dim_op
491 :
492 : CHARACTER(len=*), PARAMETER :: routineN = 'compute_berry_operator'
493 :
494 : INTEGER :: handle
495 : REAL(KIND=dp), DIMENSION(3, 6) :: vector_k
496 :
497 504 : CALL timeset(routineN, handle)
498 504 : CPASSERT(dim_op == 3 .OR. dim_op == 6)
499 504 : CPASSERT(ASSOCIATED(cell))
500 504 : vector_k = 0.0_dp
501 2016 : vector_k(:, 1) = twopi*cell%h_inv(1, :)
502 2016 : vector_k(:, 2) = twopi*cell%h_inv(2, :)
503 2016 : vector_k(:, 3) = twopi*cell%h_inv(3, :)
504 2016 : vector_k(:, 4) = twopi*(cell%h_inv(1, :) + cell%h_inv(2, :))
505 2016 : vector_k(:, 5) = twopi*(cell%h_inv(1, :) + cell%h_inv(3, :))
506 2016 : vector_k(:, 6) = twopi*(cell%h_inv(2, :) + cell%h_inv(3, :))
507 :
508 : CALL build_exp_ikr_matrix(qs_env, op_sm_set, vector_k(:, 1:dim_op), &
509 504 : force_periodic=.TRUE., cell_external=cell)
510 :
511 504 : CALL timestop(handle)
512 :
513 504 : END SUBROUTINE compute_berry_operator
514 :
515 : ! **************************************************************************************************
516 : !> \brief ...
517 : !> \param qs_loc_env ...
518 : !> \param section ...
519 : !> \param mo_array ...
520 : !> \param coeff_localized ...
521 : !> \param do_homo ...
522 : !> \param evals ...
523 : !> \param do_mixed ...
524 : ! **************************************************************************************************
525 298 : SUBROUTINE loc_write_restart(qs_loc_env, section, mo_array, coeff_localized, &
526 : do_homo, evals, do_mixed)
527 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
528 : TYPE(section_vals_type), POINTER :: section
529 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mo_array
530 : TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: coeff_localized
531 : LOGICAL, INTENT(IN) :: do_homo
532 : TYPE(cp_1d_r_p_type), DIMENSION(:), OPTIONAL, &
533 : POINTER :: evals
534 : LOGICAL, INTENT(IN), OPTIONAL :: do_mixed
535 :
536 : CHARACTER(LEN=*), PARAMETER :: routineN = 'loc_write_restart'
537 :
538 : CHARACTER(LEN=default_path_length) :: filename
539 : CHARACTER(LEN=default_string_length) :: my_middle
540 : INTEGER :: handle, ispin, max_block, nao, nloc, &
541 : nmo, output_unit, rst_unit
542 : LOGICAL :: my_do_mixed
543 : TYPE(cp_logger_type), POINTER :: logger
544 : TYPE(section_vals_type), POINTER :: print_key
545 :
546 298 : CALL timeset(routineN, handle)
547 298 : NULLIFY (logger)
548 298 : logger => cp_get_default_logger()
549 298 : output_unit = cp_logger_get_default_io_unit(logger)
550 :
551 298 : IF (qs_loc_env%do_localize) THEN
552 :
553 282 : print_key => section_vals_get_subs_vals(section, "LOC_RESTART")
554 282 : IF (BTEST(cp_print_key_should_output(logger%iter_info, &
555 : section, "LOC_RESTART"), &
556 : cp_p_file)) THEN
557 :
558 : ! Open file
559 : rst_unit = -1
560 :
561 30 : my_do_mixed = .FALSE.
562 30 : IF (PRESENT(do_mixed)) my_do_mixed = do_mixed
563 30 : IF (do_homo) THEN
564 30 : my_middle = "LOC_HOMO"
565 0 : ELSE IF (my_do_mixed) THEN
566 0 : my_middle = "LOC_MIXED"
567 : ELSE
568 0 : my_middle = "LOC_LUMO"
569 : END IF
570 :
571 : rst_unit = cp_print_key_unit_nr(logger, section, "LOC_RESTART", &
572 : extension=".wfn", file_status="REPLACE", file_action="WRITE", &
573 30 : file_form="UNFORMATTED", middle_name=TRIM(my_middle))
574 :
575 30 : IF (rst_unit > 0) filename = cp_print_key_generate_filename(logger, print_key, &
576 : middle_name=TRIM(my_middle), extension=".wfn", &
577 15 : my_local=.FALSE.)
578 :
579 30 : IF (output_unit > 0) THEN
580 : WRITE (UNIT=output_unit, FMT="(/,T2,A, A/)") &
581 15 : "LOCALIZATION| Write restart file for the localized MOS : ", &
582 30 : TRIM(filename)
583 : END IF
584 :
585 30 : IF (rst_unit > 0) THEN
586 15 : WRITE (rst_unit) qs_loc_env%localized_wfn_control%set_of_states
587 105 : WRITE (rst_unit) qs_loc_env%localized_wfn_control%lu_bound_states
588 45 : WRITE (rst_unit) qs_loc_env%localized_wfn_control%nloc_states
589 : END IF
590 :
591 70 : DO ispin = 1, SIZE(coeff_localized)
592 30 : ASSOCIATE (mo_coeff => coeff_localized(ispin))
593 40 : CALL cp_fm_get_info(mo_coeff, nrow_global=nao, ncol_global=nmo, ncol_block=max_block)
594 40 : nloc = qs_loc_env%localized_wfn_control%nloc_states(ispin)
595 40 : IF (rst_unit > 0) THEN
596 198 : WRITE (rst_unit) qs_loc_env%localized_wfn_control%loc_states(1:nloc, ispin)
597 20 : IF (do_homo .OR. my_do_mixed) THEN
598 20 : WRITE (rst_unit) nmo, &
599 20 : mo_array(ispin)%homo, &
600 20 : mo_array(ispin)%lfomo, &
601 40 : mo_array(ispin)%nelectron
602 456 : WRITE (rst_unit) mo_array(ispin)%eigenvalues(1:nmo), &
603 476 : mo_array(ispin)%occupation_numbers(1:nmo)
604 : ELSE
605 0 : WRITE (rst_unit) nmo
606 0 : WRITE (rst_unit) evals(ispin)%array(1:nmo)
607 : END IF
608 : END IF
609 :
610 80 : CALL cp_fm_write_unformatted(mo_coeff, rst_unit)
611 : END ASSOCIATE
612 :
613 : END DO
614 :
615 : CALL cp_print_key_finished_output(rst_unit, logger, section, &
616 30 : "LOC_RESTART")
617 : END IF
618 :
619 : END IF
620 :
621 298 : CALL timestop(handle)
622 :
623 298 : END SUBROUTINE loc_write_restart
624 :
625 : ! **************************************************************************************************
626 : !> \brief ...
627 : !> \param qs_loc_env ...
628 : !> \param mos ...
629 : !> \param mos_localized ...
630 : !> \param section ...
631 : !> \param section2 ...
632 : !> \param para_env ...
633 : !> \param do_homo ...
634 : !> \param restart_found ...
635 : !> \param evals ...
636 : !> \param do_mixed ...
637 : ! **************************************************************************************************
638 6 : SUBROUTINE loc_read_restart(qs_loc_env, mos, mos_localized, section, section2, para_env, &
639 : do_homo, restart_found, evals, do_mixed)
640 :
641 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
642 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
643 : TYPE(cp_fm_type), DIMENSION(:), INTENT(INOUT) :: mos_localized
644 : TYPE(section_vals_type), POINTER :: section, section2
645 : TYPE(mp_para_env_type), POINTER :: para_env
646 : LOGICAL, INTENT(IN) :: do_homo
647 : LOGICAL, INTENT(INOUT) :: restart_found
648 : TYPE(cp_1d_r_p_type), DIMENSION(:), OPTIONAL, &
649 : POINTER :: evals
650 : LOGICAL, INTENT(IN), OPTIONAL :: do_mixed
651 :
652 : CHARACTER(len=*), PARAMETER :: routineN = 'loc_read_restart'
653 :
654 : CHARACTER(LEN=25) :: fname_key
655 : CHARACTER(LEN=default_path_length) :: filename
656 : CHARACTER(LEN=default_string_length) :: my_middle
657 : INTEGER :: handle, homo_read, i, ispin, lfomo_read, max_nloc, n_rep_val, nao, &
658 : nelectron_read, nloc, nmo, nmo_read, nspin, output_unit, rst_unit
659 : LOGICAL :: file_exists, my_do_mixed
660 6 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eig_read, occ_read
661 6 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: vecbuffer
662 : TYPE(cp_logger_type), POINTER :: logger
663 : TYPE(section_vals_type), POINTER :: print_key
664 :
665 6 : CALL timeset(routineN, handle)
666 :
667 6 : logger => cp_get_default_logger()
668 :
669 6 : nspin = SIZE(mos_localized)
670 6 : nao = mos(1)%nao
671 6 : rst_unit = -1
672 :
673 : output_unit = cp_print_key_unit_nr(logger, section2, &
674 6 : "PROGRAM_RUN_INFO", extension=".Log")
675 :
676 6 : my_do_mixed = .FALSE.
677 6 : IF (PRESENT(do_mixed)) my_do_mixed = do_mixed
678 6 : IF (do_homo) THEN
679 6 : fname_key = "LOCHOMO_RESTART_FILE_NAME"
680 0 : ELSE IF (my_do_mixed) THEN
681 0 : fname_key = "LOCMIXD_RESTART_FILE_NAME"
682 : ELSE
683 0 : fname_key = "LOCLUMO_RESTART_FILE_NAME"
684 0 : IF (.NOT. PRESENT(evals)) THEN
685 0 : CPABORT("Missing argument to localize unoccupied states.")
686 : END IF
687 : END IF
688 :
689 6 : file_exists = .FALSE.
690 6 : CALL section_vals_val_get(section, fname_key, n_rep_val=n_rep_val)
691 6 : IF (n_rep_val > 0) THEN
692 0 : CALL section_vals_val_get(section, fname_key, c_val=filename)
693 : ELSE
694 :
695 6 : print_key => section_vals_get_subs_vals(section2, "LOC_RESTART")
696 6 : IF (do_homo) THEN
697 6 : my_middle = "LOC_HOMO"
698 0 : ELSE IF (my_do_mixed) THEN
699 0 : my_middle = "LOC_MIXED"
700 : ELSE
701 0 : my_middle = "LOC_LUMO"
702 : END IF
703 : filename = cp_print_key_generate_filename(logger, print_key, &
704 : middle_name=TRIM(my_middle), extension=".wfn", &
705 6 : my_local=.FALSE.)
706 : END IF
707 :
708 6 : IF (para_env%is_source()) INQUIRE (FILE=filename, exist=file_exists)
709 :
710 6 : IF (file_exists) THEN
711 2 : IF (para_env%is_source()) THEN
712 : CALL open_file(file_name=filename, &
713 : file_action="READ", &
714 : file_form="UNFORMATTED", &
715 : file_status="OLD", &
716 2 : unit_number=rst_unit)
717 :
718 2 : READ (rst_unit) qs_loc_env%localized_wfn_control%set_of_states
719 14 : READ (rst_unit) qs_loc_env%localized_wfn_control%lu_bound_states
720 6 : READ (rst_unit) qs_loc_env%localized_wfn_control%nloc_states
721 : END IF
722 : ELSE
723 4 : IF (output_unit > 0) THEN
724 : WRITE (output_unit, "(/,T10,A)") &
725 1 : "Restart file not available filename=<"//TRIM(filename)//'>'
726 : END IF
727 : END IF
728 6 : CALL para_env%bcast(file_exists)
729 :
730 6 : IF (file_exists) THEN
731 4 : restart_found = .TRUE.
732 :
733 4 : CALL para_env%bcast(qs_loc_env%localized_wfn_control%set_of_states)
734 4 : CALL para_env%bcast(qs_loc_env%localized_wfn_control%lu_bound_states)
735 4 : CALL para_env%bcast(qs_loc_env%localized_wfn_control%nloc_states)
736 :
737 12 : max_nloc = MAXVAL(qs_loc_env%localized_wfn_control%nloc_states(:))
738 :
739 12 : ALLOCATE (vecbuffer(1, nao))
740 4 : IF (ASSOCIATED(qs_loc_env%localized_wfn_control%loc_states)) THEN
741 2 : DEALLOCATE (qs_loc_env%localized_wfn_control%loc_states)
742 : END IF
743 12 : ALLOCATE (qs_loc_env%localized_wfn_control%loc_states(max_nloc, 2))
744 56 : qs_loc_env%localized_wfn_control%loc_states = 0
745 :
746 10 : DO ispin = 1, nspin
747 6 : IF (do_homo .OR. do_mixed) THEN
748 6 : nmo = mos(ispin)%nmo
749 : ELSE
750 0 : nmo = SIZE(evals(ispin)%array, 1)
751 : END IF
752 6 : IF (para_env%is_source() .AND. (nmo > 0)) THEN
753 3 : nloc = qs_loc_env%localized_wfn_control%nloc_states(ispin)
754 21 : READ (rst_unit) qs_loc_env%localized_wfn_control%loc_states(1:nloc, ispin)
755 3 : IF (do_homo .OR. do_mixed) THEN
756 3 : READ (rst_unit) nmo_read, homo_read, lfomo_read, nelectron_read
757 12 : ALLOCATE (eig_read(nmo_read), occ_read(nmo_read))
758 3 : eig_read = 0.0_dp
759 3 : occ_read = 0.0_dp
760 3 : READ (rst_unit) eig_read(1:nmo_read), occ_read(1:nmo_read)
761 : ELSE
762 0 : READ (rst_unit) nmo_read
763 0 : ALLOCATE (eig_read(nmo_read))
764 0 : eig_read = 0.0_dp
765 0 : READ (rst_unit) eig_read(1:nmo_read)
766 : END IF
767 3 : IF (nmo_read < nmo) THEN
768 : CALL cp_warn(__LOCATION__, &
769 : "The number of MOs on the restart unit is smaller than the number of "// &
770 0 : "the allocated MOs. ")
771 : END IF
772 3 : IF (nmo_read > nmo) THEN
773 : CALL cp_warn(__LOCATION__, &
774 : "The number of MOs on the restart unit is greater than the number of "// &
775 0 : "the allocated MOs. The read MO set will be truncated!")
776 : END IF
777 :
778 3 : nmo = MIN(nmo, nmo_read)
779 3 : IF (do_homo .OR. do_mixed) THEN
780 79 : mos(ispin)%eigenvalues(1:nmo) = eig_read(1:nmo)
781 79 : mos(ispin)%occupation_numbers(1:nmo) = occ_read(1:nmo)
782 3 : DEALLOCATE (eig_read, occ_read)
783 : ELSE
784 0 : evals(ispin)%array(1:nmo) = eig_read(1:nmo)
785 0 : DEALLOCATE (eig_read)
786 : END IF
787 :
788 : END IF
789 6 : IF (do_homo .OR. do_mixed) THEN
790 310 : CALL para_env%bcast(mos(ispin)%eigenvalues)
791 310 : CALL para_env%bcast(mos(ispin)%occupation_numbers)
792 : ELSE
793 0 : CALL para_env%bcast(evals(ispin)%array)
794 : END IF
795 :
796 162 : DO i = 1, nmo
797 152 : IF (para_env%is_source()) THEN
798 15236 : READ (rst_unit) vecbuffer
799 : ELSE
800 7656 : vecbuffer(1, :) = 0.0_dp
801 : END IF
802 60792 : CALL para_env%bcast(vecbuffer)
803 : CALL cp_fm_set_submatrix(mos_localized(ispin), &
804 158 : vecbuffer, 1, i, nao, 1, transpose=.TRUE.)
805 : END DO
806 : END DO
807 :
808 108 : CALL para_env%bcast(qs_loc_env%localized_wfn_control%loc_states)
809 :
810 4 : DEALLOCATE (vecbuffer)
811 :
812 : END IF
813 :
814 : ! Close restart file
815 6 : IF (para_env%is_source()) THEN
816 3 : IF (file_exists) CALL close_file(unit_number=rst_unit)
817 : END IF
818 :
819 6 : CALL timestop(handle)
820 :
821 6 : END SUBROUTINE loc_read_restart
822 :
823 : ! **************************************************************************************************
824 : !> \brief initializes everything needed for localization of the HOMOs
825 : !> \param qs_loc_env ...
826 : !> \param loc_section ...
827 : !> \param do_homo ...
828 : !> \param do_mixed ...
829 : !> \param do_xas ...
830 : !> \param nloc_xas ...
831 : !> \param spin_xas ...
832 : !> \par History
833 : !> 2009 created
834 : ! **************************************************************************************************
835 432 : SUBROUTINE qs_loc_control_init(qs_loc_env, loc_section, do_homo, do_mixed, &
836 : do_xas, nloc_xas, spin_xas)
837 :
838 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
839 : TYPE(section_vals_type), POINTER :: loc_section
840 : LOGICAL, INTENT(IN) :: do_homo
841 : LOGICAL, INTENT(IN), OPTIONAL :: do_mixed, do_xas
842 : INTEGER, INTENT(IN), OPTIONAL :: nloc_xas, spin_xas
843 :
844 : LOGICAL :: my_do_mixed
845 : TYPE(localized_wfn_control_type), POINTER :: localized_wfn_control
846 :
847 432 : NULLIFY (localized_wfn_control)
848 :
849 432 : IF (PRESENT(do_mixed)) THEN
850 2 : my_do_mixed = do_mixed
851 : ELSE
852 430 : my_do_mixed = .FALSE.
853 : END IF
854 432 : CALL localized_wfn_control_create(localized_wfn_control)
855 432 : CALL set_qs_loc_env(qs_loc_env, localized_wfn_control=localized_wfn_control)
856 432 : CALL localized_wfn_control_release(localized_wfn_control)
857 432 : CALL get_qs_loc_env(qs_loc_env, localized_wfn_control=localized_wfn_control)
858 432 : localized_wfn_control%do_homo = do_homo
859 432 : localized_wfn_control%do_mixed = my_do_mixed
860 : CALL read_loc_section(localized_wfn_control, loc_section, qs_loc_env%do_localize, &
861 432 : my_do_mixed, do_xas, nloc_xas, spin_xas)
862 :
863 432 : END SUBROUTINE qs_loc_control_init
864 :
865 : ! **************************************************************************************************
866 : !> \brief initializes everything needed for localization of the molecular orbitals
867 : !> \param qs_env ...
868 : !> \param qs_loc_env ...
869 : !> \param localize_section ...
870 : !> \param mos_localized ...
871 : !> \param do_homo ...
872 : !> \param do_mo_cubes ...
873 : !> \param mo_loc_history ...
874 : !> \param evals ...
875 : !> \param tot_zeff_corr ...
876 : !> \param do_mixed ...
877 : ! **************************************************************************************************
878 324 : SUBROUTINE qs_loc_init(qs_env, qs_loc_env, localize_section, mos_localized, &
879 : do_homo, do_mo_cubes, mo_loc_history, evals, &
880 : tot_zeff_corr, do_mixed)
881 :
882 : TYPE(qs_environment_type), POINTER :: qs_env
883 : TYPE(qs_loc_env_type), POINTER :: qs_loc_env
884 : TYPE(section_vals_type), POINTER :: localize_section
885 : TYPE(cp_fm_type), DIMENSION(:), INTENT(INOUT) :: mos_localized
886 : LOGICAL, OPTIONAL :: do_homo, do_mo_cubes
887 : TYPE(cp_fm_type), DIMENSION(:), OPTIONAL, POINTER :: mo_loc_history
888 : TYPE(cp_1d_r_p_type), DIMENSION(:), OPTIONAL, &
889 : POINTER :: evals
890 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: tot_zeff_corr
891 : LOGICAL, OPTIONAL :: do_mixed
892 :
893 : CHARACTER(len=*), PARAMETER :: routineN = 'qs_loc_init'
894 :
895 : INTEGER :: handle, homo, i, ilast_intocc, ilow, ispin, iup, n_mo(2), n_mos(2), nao, &
896 : nelectron, nextra, nmoloc(2), nocc, npocc, nspin, output_unit
897 : LOGICAL :: my_do_homo, my_do_mixed, my_do_mo_cubes, &
898 : restart_found
899 : REAL(KIND=dp) :: maxocc, my_tot_zeff_corr
900 324 : REAL(KIND=dp), DIMENSION(:), POINTER :: mo_eigenvalues, occupation
901 : TYPE(cp_fm_type), POINTER :: mo_coeff
902 : TYPE(cp_logger_type), POINTER :: logger
903 324 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: ks_rmpv, mo_derivs
904 : TYPE(dft_control_type), POINTER :: dft_control
905 : TYPE(localized_wfn_control_type), POINTER :: localized_wfn_control
906 324 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
907 : TYPE(mp_para_env_type), POINTER :: para_env
908 : TYPE(scf_control_type), POINTER :: scf_control
909 : TYPE(section_vals_type), POINTER :: loc_print_section
910 :
911 324 : CALL timeset(routineN, handle)
912 :
913 324 : NULLIFY (mos, mo_coeff, mo_eigenvalues, occupation, ks_rmpv, mo_derivs, scf_control, para_env)
914 : CALL get_qs_env(qs_env, &
915 : mos=mos, &
916 : matrix_ks=ks_rmpv, &
917 : mo_derivs=mo_derivs, &
918 : scf_control=scf_control, &
919 : dft_control=dft_control, &
920 324 : para_env=para_env)
921 :
922 324 : loc_print_section => section_vals_get_subs_vals(localize_section, "PRINT")
923 :
924 324 : logger => cp_get_default_logger()
925 324 : output_unit = cp_logger_get_default_io_unit(logger)
926 :
927 324 : nspin = SIZE(mos)
928 324 : IF (PRESENT(do_homo)) THEN
929 324 : my_do_homo = do_homo
930 : ELSE
931 0 : my_do_homo = .TRUE.
932 : END IF
933 324 : IF (PRESENT(do_mo_cubes)) THEN
934 134 : my_do_mo_cubes = do_mo_cubes
935 : ELSE
936 : my_do_mo_cubes = .FALSE.
937 : END IF
938 324 : IF (PRESENT(do_mixed)) THEN
939 2 : my_do_mixed = do_mixed
940 : ELSE
941 322 : my_do_mixed = .FALSE.
942 : END IF
943 324 : IF (PRESENT(tot_zeff_corr)) THEN
944 2 : my_tot_zeff_corr = tot_zeff_corr
945 : ELSE
946 322 : my_tot_zeff_corr = 0.0_dp
947 : END IF
948 324 : restart_found = .FALSE.
949 :
950 324 : IF (qs_loc_env%do_localize) THEN
951 : ! Some setup for MOs to be localized
952 308 : CALL get_qs_loc_env(qs_loc_env, localized_wfn_control=localized_wfn_control)
953 308 : IF (localized_wfn_control%loc_restart) THEN
954 6 : IF (localized_wfn_control%nextra > 0) THEN
955 : ! currently only the occupied guess is read
956 0 : my_do_homo = .FALSE.
957 : END IF
958 : CALL loc_read_restart(qs_loc_env, mos, mos_localized, localize_section, &
959 : loc_print_section, para_env, my_do_homo, restart_found, evals=evals, &
960 6 : do_mixed=my_do_mixed)
961 9 : IF (output_unit > 0) WRITE (output_unit, "(/,T2,A,A)") "LOCALIZATION| ", &
962 6 : " The orbitals to be localized are read from localization restart file."
963 18 : nmoloc = localized_wfn_control%nloc_states
964 18 : localized_wfn_control%nguess = nmoloc
965 6 : IF (localized_wfn_control%nextra > 0) THEN
966 : ! reset different variables in localized_wfn_control:
967 : ! lu_bound_states, nloc_states, loc_states
968 0 : localized_wfn_control%loc_restart = restart_found
969 0 : localized_wfn_control%set_of_states = state_loc_mixed
970 0 : DO ispin = 1, nspin
971 : CALL get_mo_set(mos(ispin), homo=homo, occupation_numbers=occupation, &
972 0 : maxocc=maxocc)
973 0 : nextra = localized_wfn_control%nextra
974 0 : nocc = homo
975 0 : DO i = nocc, 1, -1
976 0 : IF (maxocc - occupation(i) < localized_wfn_control%eps_occ) THEN
977 0 : ilast_intocc = i
978 0 : EXIT
979 : END IF
980 : END DO
981 0 : nocc = ilast_intocc
982 0 : npocc = homo - nocc
983 0 : nmoloc(ispin) = nocc + nextra
984 0 : localized_wfn_control%lu_bound_states(1, ispin) = 1
985 0 : localized_wfn_control%lu_bound_states(2, ispin) = nmoloc(ispin)
986 0 : localized_wfn_control%nloc_states(ispin) = nmoloc(ispin)
987 : END DO
988 0 : my_do_homo = .FALSE.
989 : END IF
990 : END IF
991 308 : IF (.NOT. restart_found) THEN
992 304 : nmoloc = 0
993 726 : DO ispin = 1, nspin
994 : CALL get_mo_set(mos(ispin), nmo=n_mo(ispin), nelectron=nelectron, homo=homo, nao=nao, &
995 : mo_coeff=mo_coeff, eigenvalues=mo_eigenvalues, occupation_numbers=occupation, &
996 422 : maxocc=maxocc)
997 : ! Get eigenstates (only needed if not already calculated before)
998 : IF ((.NOT. my_do_mo_cubes) &
999 : .AND. my_do_homo .AND. ASSOCIATED(qs_env%scf_env) &
1000 422 : .AND. qs_env%scf_env%method == ot_method_nr .AND. (.NOT. dft_control%restricted)) THEN
1001 30 : CALL make_mo_eig(mos, nspin, ks_rmpv, scf_control, mo_derivs)
1002 : END IF
1003 1148 : IF (localized_wfn_control%set_of_states == state_loc_all .AND. my_do_homo) THEN
1004 382 : nmoloc(ispin) = NINT(nelectron/occupation(1))
1005 382 : IF (n_mo(ispin) > homo) THEN
1006 28 : DO i = nmoloc(ispin), 1, -1
1007 28 : IF (occupation(1) - occupation(i) < localized_wfn_control%eps_occ) THEN
1008 14 : ilast_intocc = i
1009 14 : EXIT
1010 : END IF
1011 : END DO
1012 : ELSE
1013 368 : ilast_intocc = nmoloc(ispin)
1014 : END IF
1015 382 : nmoloc(ispin) = ilast_intocc
1016 382 : localized_wfn_control%lu_bound_states(1, ispin) = 1
1017 382 : localized_wfn_control%lu_bound_states(2, ispin) = ilast_intocc
1018 382 : IF (nmoloc(ispin) /= n_mo(ispin)) THEN
1019 14 : IF (output_unit > 0) THEN
1020 : WRITE (output_unit, "(/,T2,A,I4,A,I6,A,/,T15,A,F12.6,A,F12.6,A)") &
1021 7 : "LOCALIZATION| Spin ", ispin, " The first ", &
1022 7 : ilast_intocc, " occupied orbitals are localized,", " with energies from ", &
1023 14 : mo_eigenvalues(1), " to ", mo_eigenvalues(ilast_intocc), " [a.u.]."
1024 : END IF
1025 : END IF
1026 40 : ELSE IF (localized_wfn_control%set_of_states == energy_loc_range .AND. my_do_homo) THEN
1027 12 : ilow = 0
1028 12 : iup = 0
1029 20 : DO i = 1, n_mo(ispin)
1030 20 : IF (mo_eigenvalues(i) >= localized_wfn_control%lu_ene_bound(1)) THEN
1031 : ilow = i
1032 : EXIT
1033 : END IF
1034 : END DO
1035 306 : DO i = n_mo(ispin), 1, -1
1036 306 : IF (mo_eigenvalues(i) <= localized_wfn_control%lu_ene_bound(2)) THEN
1037 : iup = i
1038 : EXIT
1039 : END IF
1040 : END DO
1041 12 : localized_wfn_control%lu_bound_states(1, ispin) = ilow
1042 12 : localized_wfn_control%lu_bound_states(2, ispin) = iup
1043 12 : localized_wfn_control%nloc_states(ispin) = iup - ilow + 1
1044 12 : nmoloc(ispin) = localized_wfn_control%nloc_states(ispin)
1045 12 : IF (occupation(ilow) - occupation(iup) > localized_wfn_control%eps_occ) THEN
1046 : CALL cp_abort(__LOCATION__, &
1047 : "The selected energy range includes orbitals with different occupation number. "// &
1048 0 : "The localization procedure cannot be applied.")
1049 : END IF
1050 18 : IF (output_unit > 0) WRITE (output_unit, "(/,T2,A,I4,A,I6,A)") "LOCALIZATION| Spin ", ispin, " : ", &
1051 12 : nmoloc(ispin), " orbitals in the selected energy range are localized."
1052 28 : ELSE IF (localized_wfn_control%set_of_states == state_loc_all .AND. (.NOT. my_do_homo)) THEN
1053 0 : nmoloc(ispin) = n_mo(ispin) - homo
1054 0 : localized_wfn_control%lu_bound_states(1, ispin) = homo + 1
1055 0 : localized_wfn_control%lu_bound_states(2, ispin) = n_mo(ispin)
1056 0 : IF (output_unit > 0) THEN
1057 : WRITE (output_unit, "(/,T2,A,I4,A,I6,A,/,T15,A,F12.6,A,F12.6,A)") &
1058 0 : "LOCALIZATION| Spin ", ispin, " The first ", &
1059 0 : nmoloc(ispin), " virtual orbitals are localized,", " with energies from ", &
1060 0 : mo_eigenvalues(homo + 1), " to ", mo_eigenvalues(n_mo(ispin)), " [a.u.]."
1061 : END IF
1062 28 : ELSE IF (localized_wfn_control%set_of_states == state_loc_mixed) THEN
1063 2 : nextra = localized_wfn_control%nextra
1064 2 : nocc = homo
1065 6 : DO i = nocc, 1, -1
1066 6 : IF (maxocc - occupation(i) < localized_wfn_control%eps_occ) THEN
1067 2 : ilast_intocc = i
1068 2 : EXIT
1069 : END IF
1070 : END DO
1071 2 : nocc = ilast_intocc
1072 2 : npocc = homo - nocc
1073 2 : nmoloc(ispin) = nocc + nextra
1074 2 : localized_wfn_control%lu_bound_states(1, ispin) = 1
1075 2 : localized_wfn_control%lu_bound_states(2, ispin) = nmoloc(ispin)
1076 2 : IF (output_unit > 0) THEN
1077 : WRITE (output_unit, "(/,T2,A,I4,A,I6,A,/,T15,A,I6,/,T15,A,I6,/,T15,A,I6,/,T15,A,F12.6,A)") &
1078 1 : "LOCALIZATION| Spin ", ispin, " The first ", &
1079 1 : nmoloc(ispin), " orbitals are localized.", &
1080 1 : "Number of fully occupied MOs: ", nocc, &
1081 1 : "Number of partially occupied MOs: ", npocc, &
1082 1 : "Number of extra degrees of freedom: ", nextra, &
1083 2 : "Excess charge: ", my_tot_zeff_corr, " electrons"
1084 : END IF
1085 : ELSE
1086 26 : nmoloc(ispin) = MIN(localized_wfn_control%nloc_states(1), n_mo(ispin))
1087 33 : IF (output_unit > 0 .AND. my_do_homo) WRITE (output_unit, "(/,T2,A,I4,A,I6,A)") "LOCALIZATION| Spin ", ispin, &
1088 14 : " : ", nmoloc(ispin), " occupied orbitals are localized, as given in the input list."
1089 19 : IF (output_unit > 0 .AND. (.NOT. my_do_homo)) WRITE (output_unit, "(/,T2,A,I4,A,I6,A)") "LOCALIZATION| Spin ", &
1090 12 : ispin, " : ", nmoloc(ispin), " unoccupied orbitals are localized, as given in the input list."
1091 26 : IF (n_mo(ispin) > homo .AND. my_do_homo) THEN
1092 8 : ilow = localized_wfn_control%loc_states(1, ispin)
1093 56 : DO i = 2, nmoloc(ispin)
1094 48 : iup = localized_wfn_control%loc_states(i, ispin)
1095 56 : IF (ABS(occupation(ilow) - occupation(iup)) > localized_wfn_control%eps_occ) THEN
1096 : ! write warning
1097 : CALL cp_warn(__LOCATION__, &
1098 : "User requested the calculation of localized wavefunction from a subset of MOs, "// &
1099 : "including MOs with different occupations. Check the selected subset, "// &
1100 : "the electronic density is not invariant with "// &
1101 0 : "respect to rotations among orbitals with different occupation numbers!")
1102 : END IF
1103 : END DO
1104 : END IF
1105 : END IF
1106 : END DO ! ispin
1107 912 : n_mos(:) = nao - n_mo(:)
1108 304 : IF (my_do_homo .OR. my_do_mixed) n_mos = n_mo
1109 304 : CALL set_loc_wfn_lists(localized_wfn_control, nmoloc, n_mos, nspin)
1110 : END IF
1111 308 : CALL set_loc_centers(localized_wfn_control, nmoloc, nspin)
1112 308 : IF (my_do_homo .OR. my_do_mixed) THEN
1113 : CALL qs_loc_env_init(qs_loc_env, localized_wfn_control, qs_env, &
1114 302 : loc_coeff=mos_localized, mo_loc_history=mo_loc_history)
1115 : END IF
1116 : ELSE
1117 : ! Let's inform in case the section is not present in the input
1118 : CALL cp_warn(__LOCATION__, &
1119 : "User requested the calculation of the localized wavefunction but the section "// &
1120 16 : "LOCALIZE was not specified. Localization will not be performed!")
1121 : END IF
1122 :
1123 324 : CALL timestop(handle)
1124 :
1125 324 : END SUBROUTINE qs_loc_init
1126 :
1127 : ! **************************************************************************************************
1128 : !> \brief read the controlparameter from input, using the new input scheme
1129 : !> \param localized_wfn_control ...
1130 : !> \param loc_section ...
1131 : !> \param localize ...
1132 : !> \param do_mixed ...
1133 : !> \param do_xas ...
1134 : !> \param nloc_xas ...
1135 : !> \param spin_channel_xas ...
1136 : !> \par History
1137 : !> 05.2005 created [MI]
1138 : ! **************************************************************************************************
1139 864 : SUBROUTINE read_loc_section(localized_wfn_control, loc_section, &
1140 : localize, do_mixed, do_xas, nloc_xas, spin_channel_xas)
1141 :
1142 : TYPE(localized_wfn_control_type), POINTER :: localized_wfn_control
1143 : TYPE(section_vals_type), POINTER :: loc_section
1144 : LOGICAL, INTENT(OUT) :: localize
1145 : LOGICAL, INTENT(IN), OPTIONAL :: do_mixed, do_xas
1146 : INTEGER, INTENT(IN), OPTIONAL :: nloc_xas, spin_channel_xas
1147 :
1148 : INTEGER :: i, ind, ir, n_list, n_rep, n_state, &
1149 : nextra, nline, other_spin, &
1150 : output_unit, spin_xas
1151 432 : INTEGER, DIMENSION(:), POINTER :: list, loc_list
1152 : LOGICAL :: my_do_mixed, my_do_xas
1153 432 : REAL(dp), POINTER :: ene(:)
1154 : TYPE(cp_logger_type), POINTER :: logger
1155 : TYPE(section_vals_type), POINTER :: loc_print_section
1156 :
1157 432 : my_do_xas = .FALSE.
1158 432 : spin_xas = 1
1159 432 : IF (PRESENT(do_xas)) THEN
1160 108 : my_do_xas = do_xas
1161 108 : CPASSERT(PRESENT(nloc_xas))
1162 : END IF
1163 432 : IF (PRESENT(spin_channel_xas)) spin_xas = spin_channel_xas
1164 432 : my_do_mixed = .FALSE.
1165 432 : IF (PRESENT(do_mixed)) THEN
1166 432 : my_do_mixed = do_mixed
1167 : END IF
1168 432 : CPASSERT(ASSOCIATED(loc_section))
1169 432 : NULLIFY (logger)
1170 432 : logger => cp_get_default_logger()
1171 :
1172 432 : CALL section_vals_val_get(loc_section, "_SECTION_PARAMETERS_", l_val=localize)
1173 432 : IF (localize) THEN
1174 384 : loc_print_section => section_vals_get_subs_vals(loc_section, "PRINT")
1175 384 : NULLIFY (list)
1176 384 : NULLIFY (loc_list)
1177 2688 : localized_wfn_control%lu_bound_states = 0
1178 1152 : localized_wfn_control%lu_ene_bound = 0.0_dp
1179 1152 : localized_wfn_control%nloc_states = 0
1180 384 : localized_wfn_control%set_of_states = 0
1181 384 : localized_wfn_control%nextra = 0
1182 384 : n_state = 0
1183 :
1184 : CALL section_vals_val_get(loc_section, "MAX_ITER", &
1185 384 : i_val=localized_wfn_control%max_iter)
1186 : CALL section_vals_val_get(loc_section, "MAX_CRAZY_ANGLE", &
1187 384 : r_val=localized_wfn_control%max_crazy_angle)
1188 : CALL section_vals_val_get(loc_section, "CRAZY_SCALE", &
1189 384 : r_val=localized_wfn_control%crazy_scale)
1190 : CALL section_vals_val_get(loc_section, "EPS_OCCUPATION", &
1191 384 : r_val=localized_wfn_control%eps_occ)
1192 : CALL section_vals_val_get(loc_section, "CRAZY_USE_DIAG", &
1193 384 : l_val=localized_wfn_control%crazy_use_diag)
1194 : CALL section_vals_val_get(loc_section, "OUT_ITER_EACH", &
1195 384 : i_val=localized_wfn_control%out_each)
1196 : CALL section_vals_val_get(loc_section, "EPS_LOCALIZATION", &
1197 384 : r_val=localized_wfn_control%eps_localization)
1198 : CALL section_vals_val_get(loc_section, "MIN_OR_MAX", &
1199 384 : i_val=localized_wfn_control%min_or_max)
1200 : CALL section_vals_val_get(loc_section, "JACOBI_FALLBACK", &
1201 384 : l_val=localized_wfn_control%jacobi_fallback)
1202 : CALL section_vals_val_get(loc_section, "JACOBI_REFINEMENT", &
1203 384 : l_val=localized_wfn_control%jacobi_refinement)
1204 : CALL section_vals_val_get(loc_section, "METHOD", &
1205 384 : i_val=localized_wfn_control%localization_method)
1206 : CALL section_vals_val_get(loc_section, "OPERATOR", &
1207 384 : i_val=localized_wfn_control%operator_type)
1208 : CALL section_vals_val_get(loc_section, "RESTART", &
1209 384 : l_val=localized_wfn_control%loc_restart)
1210 : CALL section_vals_val_get(loc_section, "USE_HISTORY", &
1211 384 : l_val=localized_wfn_control%use_history)
1212 : CALL section_vals_val_get(loc_section, "NEXTRA", &
1213 384 : i_val=localized_wfn_control%nextra)
1214 : CALL section_vals_val_get(loc_section, "CPO_GUESS", &
1215 384 : i_val=localized_wfn_control%coeff_po_guess)
1216 : CALL section_vals_val_get(loc_section, "CPO_GUESS_SPACE", &
1217 384 : i_val=localized_wfn_control%coeff_po_guess_mo_space)
1218 : CALL section_vals_val_get(loc_section, "CG_PO", &
1219 384 : l_val=localized_wfn_control%do_cg_po)
1220 :
1221 384 : IF (localized_wfn_control%do_homo) THEN
1222 : ! List of States HOMO
1223 376 : CALL section_vals_val_get(loc_section, "LIST", n_rep_val=n_rep)
1224 376 : IF (n_rep > 0) THEN
1225 14 : n_list = 0
1226 28 : DO ir = 1, n_rep
1227 14 : NULLIFY (list)
1228 14 : CALL section_vals_val_get(loc_section, "LIST", i_rep_val=ir, i_vals=list)
1229 28 : IF (ASSOCIATED(list)) THEN
1230 14 : CALL reallocate(loc_list, 1, n_list + SIZE(list))
1231 90 : DO i = 1, SIZE(list)
1232 90 : loc_list(n_list + i) = list(i)
1233 : END DO ! i
1234 14 : n_list = n_list + SIZE(list)
1235 : END IF
1236 : END DO ! ir
1237 14 : IF (n_list /= 0) THEN
1238 14 : localized_wfn_control%set_of_states = state_loc_list
1239 42 : ALLOCATE (localized_wfn_control%loc_states(n_list, 2))
1240 194 : localized_wfn_control%loc_states = 0
1241 180 : localized_wfn_control%loc_states(:, 1) = loc_list(:)
1242 180 : localized_wfn_control%loc_states(:, 2) = loc_list(:)
1243 14 : localized_wfn_control%nloc_states(1) = n_list
1244 14 : localized_wfn_control%nloc_states(2) = n_list
1245 14 : IF (my_do_xas) THEN
1246 4 : other_spin = 2
1247 4 : IF (spin_xas == 2) other_spin = 1
1248 4 : localized_wfn_control%nloc_states(other_spin) = 0
1249 22 : localized_wfn_control%loc_states(:, other_spin) = 0
1250 : END IF
1251 14 : DEALLOCATE (loc_list)
1252 : END IF
1253 : END IF
1254 :
1255 : ELSE
1256 : ! List of States LUMO
1257 8 : CALL section_vals_val_get(loc_section, "LIST_UNOCCUPIED", n_rep_val=n_rep)
1258 8 : IF (n_rep > 0) THEN
1259 6 : n_list = 0
1260 12 : DO ir = 1, n_rep
1261 6 : NULLIFY (list)
1262 6 : CALL section_vals_val_get(loc_section, "LIST_UNOCCUPIED", i_rep_val=ir, i_vals=list)
1263 12 : IF (ASSOCIATED(list)) THEN
1264 6 : CALL reallocate(loc_list, 1, n_list + SIZE(list))
1265 46 : DO i = 1, SIZE(list)
1266 46 : loc_list(n_list + i) = list(i)
1267 : END DO ! i
1268 6 : n_list = n_list + SIZE(list)
1269 : END IF
1270 : END DO ! ir
1271 6 : IF (n_list /= 0) THEN
1272 6 : localized_wfn_control%set_of_states = state_loc_list
1273 18 : ALLOCATE (localized_wfn_control%loc_states(n_list, 2))
1274 98 : localized_wfn_control%loc_states = 0
1275 92 : localized_wfn_control%loc_states(:, 1) = loc_list(:)
1276 92 : localized_wfn_control%loc_states(:, 2) = loc_list(:)
1277 6 : localized_wfn_control%nloc_states(1) = n_list
1278 6 : DEALLOCATE (loc_list)
1279 : END IF
1280 : END IF
1281 : END IF
1282 :
1283 384 : IF (localized_wfn_control%set_of_states == 0) THEN
1284 364 : CALL section_vals_val_get(loc_section, "ENERGY_RANGE", r_vals=ene)
1285 364 : IF (ene(1) /= ene(2)) THEN
1286 10 : localized_wfn_control%set_of_states = energy_loc_range
1287 10 : localized_wfn_control%lu_ene_bound(1) = ene(1)
1288 10 : localized_wfn_control%lu_ene_bound(2) = ene(2)
1289 : END IF
1290 : END IF
1291 :
1292 : ! All States or XAS specific states
1293 384 : IF (localized_wfn_control%set_of_states == 0) THEN
1294 354 : IF (my_do_xas) THEN
1295 72 : localized_wfn_control%set_of_states = state_loc_range
1296 216 : localized_wfn_control%nloc_states(:) = 0
1297 216 : localized_wfn_control%lu_bound_states(1, :) = 0
1298 216 : localized_wfn_control%lu_bound_states(2, :) = 0
1299 72 : localized_wfn_control%nloc_states(spin_xas) = nloc_xas
1300 72 : localized_wfn_control%lu_bound_states(1, spin_xas) = 1
1301 72 : localized_wfn_control%lu_bound_states(2, spin_xas) = nloc_xas
1302 282 : ELSE IF (my_do_mixed) THEN
1303 2 : localized_wfn_control%set_of_states = state_loc_mixed
1304 2 : nextra = localized_wfn_control%nextra
1305 : ELSE
1306 280 : localized_wfn_control%set_of_states = state_loc_all
1307 : END IF
1308 : END IF
1309 :
1310 : localized_wfn_control%print_centers = &
1311 : BTEST(cp_print_key_should_output(logger%iter_info, loc_print_section, &
1312 384 : "WANNIER_CENTERS"), cp_p_file)
1313 : localized_wfn_control%print_spreads = &
1314 : BTEST(cp_print_key_should_output(logger%iter_info, loc_print_section, &
1315 384 : "WANNIER_SPREADS"), cp_p_file)
1316 : localized_wfn_control%print_cubes = &
1317 : BTEST(cp_print_key_should_output(logger%iter_info, loc_print_section, &
1318 384 : "WANNIER_CUBES"), cp_p_file)
1319 :
1320 : output_unit = cp_print_key_unit_nr(logger, loc_print_section, "PROGRAM_RUN_INFO", &
1321 384 : extension=".Log")
1322 :
1323 384 : IF (output_unit > 0) THEN
1324 : WRITE (UNIT=output_unit, FMT="(/,T2,A)") &
1325 192 : "LOCALIZE| The spread relative to a set of orbitals is computed"
1326 :
1327 332 : SELECT CASE (localized_wfn_control%set_of_states)
1328 : CASE (state_loc_all)
1329 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1330 140 : "LOCALIZE| Orbitals to be localized: All orbitals"
1331 : WRITE (UNIT=output_unit, FMT="(T2,A,/,T12,A,F16.8)") &
1332 140 : "LOCALIZE| If fractional occupation, fully occupied MOs are those ", &
1333 280 : "within occupation tolerance of ", localized_wfn_control%eps_occ
1334 : CASE (state_loc_range)
1335 : WRITE (UNIT=output_unit, FMT="(T2,A,T65,I8,A,I8)") &
1336 36 : "LOCALIZE| Orbitals to be localized: Those with index between ", &
1337 36 : localized_wfn_control%lu_bound_states(1, spin_xas), " and ", &
1338 72 : localized_wfn_control%lu_bound_states(2, spin_xas)
1339 : CASE (state_loc_list)
1340 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1341 10 : "LOCALIZE| Orbitals to be localized: Those with index in the following list"
1342 10 : nline = localized_wfn_control%nloc_states(1)/10 + 1
1343 10 : ind = 0
1344 21 : DO i = 1, nline
1345 21 : IF (ind + 10 < localized_wfn_control%nloc_states(1)) THEN
1346 11 : WRITE (UNIT=output_unit, FMT="(T8,10I7)") localized_wfn_control%loc_states(ind + 1:ind + 10, 1)
1347 1 : ind = ind + 10
1348 : ELSE
1349 : WRITE (UNIT=output_unit, FMT="(T8,10I7)") &
1350 58 : localized_wfn_control%loc_states(ind + 1:localized_wfn_control%nloc_states(1), 1)
1351 10 : ind = localized_wfn_control%nloc_states(1)
1352 : END IF
1353 : END DO
1354 : CASE (energy_loc_range)
1355 : WRITE (UNIT=output_unit, FMT="(T2,A,T65,/,f16.6,A,f16.6,A)") &
1356 5 : "LOCALIZE| Orbitals to be localized: Those with energy in the range between ", &
1357 10 : localized_wfn_control%lu_ene_bound(1), " and ", localized_wfn_control%lu_ene_bound(2), " a.u."
1358 : CASE (state_loc_mixed)
1359 : WRITE (UNIT=output_unit, FMT="(T2,A,I4,A)") &
1360 1 : "LOCALIZE| Orbitals to be localized: Occupied orbitals + ", nextra, " orbitals"
1361 : CASE DEFAULT
1362 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1363 192 : "LOCALIZE| Orbitals to be localized: None "
1364 : END SELECT
1365 :
1366 381 : SELECT CASE (localized_wfn_control%operator_type)
1367 : CASE (op_loc_berry)
1368 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1369 189 : "LOCALIZE| Spread defined by the Berry phase operator "
1370 : CASE (op_loc_boys)
1371 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1372 0 : "LOCALIZE| Spread defined by the Boys phase operator "
1373 : CASE DEFAULT
1374 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1375 192 : "LOCALIZE| Spread defined by the Pipek phase operator "
1376 : END SELECT
1377 :
1378 328 : SELECT CASE (localized_wfn_control%localization_method)
1379 : CASE (do_loc_jacobi)
1380 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1381 136 : "LOCALIZE| Optimal unitary transformation generated by Jacobi algorithm"
1382 : CASE (do_loc_crazy)
1383 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1384 38 : "LOCALIZE| Optimal unitary transformation generated by Crazy angle algorithm"
1385 : WRITE (UNIT=output_unit, FMT="(T2,A,F16.8)") &
1386 38 : "LOCALIZE| maximum angle: ", localized_wfn_control%max_crazy_angle
1387 : WRITE (UNIT=output_unit, FMT="(T2,A,F16.8)") &
1388 38 : "LOCALIZE| scaling: ", localized_wfn_control%crazy_scale
1389 : WRITE (UNIT=output_unit, FMT="(T2,A,L1)") &
1390 38 : "LOCALIZE| use diag:", localized_wfn_control%crazy_use_diag
1391 : CASE (do_loc_gapo)
1392 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1393 1 : "LOCALIZE| Optimal unitary transformation generated by gradient ascent algorithm "
1394 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1395 1 : "LOCALIZE| for partially occupied wannier functions"
1396 : CASE (do_loc_direct)
1397 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1398 1 : "LOCALIZE| Optimal unitary transformation generated by direct algorithm"
1399 : CASE (do_loc_l1_norm_sd)
1400 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1401 9 : "LOCALIZE| Optimal unitary transformation generated by "
1402 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1403 9 : "LOCALIZE| steepest descent algorithm applied on an approximate l1 norm"
1404 : CASE (do_loc_none)
1405 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1406 0 : "LOCALIZE| No unitary transformation is applied"
1407 : CASE (do_loc_scdm)
1408 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1409 192 : "LOCALIZE| Pivoted QR decomposition is used to transform coefficients"
1410 : END SELECT
1411 :
1412 : END IF ! process has output_unit
1413 :
1414 384 : CALL cp_print_key_finished_output(output_unit, logger, loc_print_section, "PROGRAM_RUN_INFO")
1415 :
1416 : ELSE
1417 48 : localized_wfn_control%localization_method = do_loc_none
1418 48 : localized_wfn_control%localization_method = state_loc_none
1419 48 : localized_wfn_control%print_centers = .FALSE.
1420 48 : localized_wfn_control%print_spreads = .FALSE.
1421 48 : localized_wfn_control%print_cubes = .FALSE.
1422 : END IF
1423 :
1424 432 : END SUBROUTINE read_loc_section
1425 :
1426 : ! **************************************************************************************************
1427 : !> \brief create the center and spread array and the file names for the output
1428 : !> \param localized_wfn_control ...
1429 : !> \param nmoloc ...
1430 : !> \param nspins ...
1431 : !> \par History
1432 : !> 04.2005 created [MI]
1433 : ! **************************************************************************************************
1434 484 : SUBROUTINE set_loc_centers(localized_wfn_control, nmoloc, nspins)
1435 :
1436 : TYPE(localized_wfn_control_type) :: localized_wfn_control
1437 : INTEGER, DIMENSION(2), INTENT(IN) :: nmoloc
1438 : INTEGER, INTENT(IN) :: nspins
1439 :
1440 : INTEGER :: ispin
1441 :
1442 1144 : DO ispin = 1, nspins
1443 1938 : ALLOCATE (localized_wfn_control%centers_set(ispin)%array(6, nmoloc(ispin)))
1444 32182 : localized_wfn_control%centers_set(ispin)%array = 0.0_dp
1445 : END DO
1446 :
1447 484 : END SUBROUTINE set_loc_centers
1448 :
1449 : ! **************************************************************************************************
1450 : !> \brief create the lists of mos that are taken into account
1451 : !> \param localized_wfn_control ...
1452 : !> \param nmoloc ...
1453 : !> \param nmo ...
1454 : !> \param nspins ...
1455 : !> \param my_spin ...
1456 : !> \par History
1457 : !> 04.2005 created [MI]
1458 : ! **************************************************************************************************
1459 346 : SUBROUTINE set_loc_wfn_lists(localized_wfn_control, nmoloc, nmo, nspins, my_spin)
1460 :
1461 : TYPE(localized_wfn_control_type) :: localized_wfn_control
1462 : INTEGER, DIMENSION(2), INTENT(IN) :: nmoloc, nmo
1463 : INTEGER, INTENT(IN) :: nspins
1464 : INTEGER, INTENT(IN), OPTIONAL :: my_spin
1465 :
1466 : CHARACTER(len=*), PARAMETER :: routineN = 'set_loc_wfn_lists'
1467 :
1468 : INTEGER :: i, ispin, max_iloc, max_nmoloc, state
1469 :
1470 346 : CALL timeset(routineN, state)
1471 :
1472 1038 : localized_wfn_control%nloc_states(1:2) = nmoloc(1:2)
1473 346 : max_nmoloc = MAX(nmoloc(1), nmoloc(2))
1474 :
1475 364 : SELECT CASE (localized_wfn_control%set_of_states)
1476 : CASE (state_loc_list)
1477 : ! List
1478 18 : CPASSERT(ASSOCIATED(localized_wfn_control%loc_states))
1479 52 : DO ispin = 1, nspins
1480 34 : localized_wfn_control%lu_bound_states(1, ispin) = 1
1481 34 : localized_wfn_control%lu_bound_states(2, ispin) = nmoloc(ispin)
1482 52 : IF (nmoloc(ispin) < 1) THEN
1483 4 : localized_wfn_control%lu_bound_states(1, ispin) = 0
1484 22 : localized_wfn_control%loc_states(:, ispin) = 0
1485 : END IF
1486 : END DO
1487 : CASE (state_loc_range)
1488 : ! Range
1489 114 : ALLOCATE (localized_wfn_control%loc_states(max_nmoloc, 2))
1490 446 : localized_wfn_control%loc_states = 0
1491 114 : DO ispin = 1, nspins
1492 : localized_wfn_control%lu_bound_states(1, ispin) = &
1493 76 : localized_wfn_control%lu_bound_states(1, my_spin)
1494 : localized_wfn_control%lu_bound_states(2, ispin) = &
1495 76 : localized_wfn_control%lu_bound_states(1, my_spin) + nmoloc(ispin) - 1
1496 76 : max_iloc = localized_wfn_control%lu_bound_states(2, ispin)
1497 242 : DO i = 1, nmoloc(ispin)
1498 242 : localized_wfn_control%loc_states(i, ispin) = localized_wfn_control%lu_bound_states(1, ispin) + i - 1
1499 : END DO
1500 76 : CPASSERT(max_iloc <= nmo(ispin))
1501 38 : MARK_USED(nmo)
1502 : END DO
1503 : CASE (energy_loc_range)
1504 : ! Energy
1505 30 : ALLOCATE (localized_wfn_control%loc_states(max_nmoloc, 2))
1506 202 : localized_wfn_control%loc_states = 0
1507 22 : DO ispin = 1, nspins
1508 128 : DO i = 1, nmoloc(ispin)
1509 118 : localized_wfn_control%loc_states(i, ispin) = localized_wfn_control%lu_bound_states(1, ispin) + i - 1
1510 : END DO
1511 : END DO
1512 : CASE (state_loc_all)
1513 : ! All
1514 834 : ALLOCATE (localized_wfn_control%loc_states(max_nmoloc, 2))
1515 5678 : localized_wfn_control%loc_states = 0
1516 :
1517 278 : IF (localized_wfn_control%lu_bound_states(1, 1) == 1) THEN
1518 660 : DO ispin = 1, nspins
1519 382 : localized_wfn_control%lu_bound_states(1, ispin) = 1
1520 382 : localized_wfn_control%lu_bound_states(2, ispin) = nmoloc(ispin)
1521 382 : IF (nmoloc(ispin) < 1) localized_wfn_control%lu_bound_states(1, ispin) = 0
1522 3808 : DO i = 1, nmoloc(ispin)
1523 3530 : localized_wfn_control%loc_states(i, ispin) = i
1524 : END DO
1525 : END DO
1526 : ELSE
1527 0 : DO ispin = 1, nspins
1528 0 : IF (nmoloc(ispin) < 1) localized_wfn_control%lu_bound_states(1, ispin) = 0
1529 0 : DO i = 1, nmoloc(ispin)
1530 : localized_wfn_control%loc_states(i, ispin) = &
1531 0 : localized_wfn_control%lu_bound_states(1, ispin) + i - 1
1532 : END DO
1533 : END DO
1534 : END IF
1535 : CASE (state_loc_mixed)
1536 : ! Mixed
1537 6 : ALLOCATE (localized_wfn_control%loc_states(max_nmoloc, 2))
1538 178 : localized_wfn_control%loc_states = 0
1539 350 : DO ispin = 1, nspins
1540 90 : DO i = 1, nmoloc(ispin)
1541 88 : localized_wfn_control%loc_states(i, ispin) = i
1542 : END DO
1543 : END DO
1544 : END SELECT
1545 :
1546 346 : CALL timestop(state)
1547 :
1548 346 : END SUBROUTINE set_loc_wfn_lists
1549 :
1550 : END MODULE qs_loc_utils
|
