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 : MODULE optbas_fenv_manipulation
8 : USE atomic_kind_types, ONLY: atomic_kind_type,&
9 : get_atomic_kind
10 : USE basis_set_container_types, ONLY: get_basis_from_container
11 : USE basis_set_types, ONLY: gto_basis_set_type,&
12 : init_orb_basis_set
13 : USE cp_blacs_env, ONLY: cp_blacs_env_type
14 : USE cp_control_types, ONLY: dft_control_type
15 : USE cp_dbcsr_api, ONLY: dbcsr_get_info,&
16 : dbcsr_p_type,&
17 : dbcsr_type
18 : USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm
19 : USE cp_fm_basic_linalg, ONLY: cp_fm_uplo_to_full
20 : USE cp_fm_cholesky, ONLY: cp_fm_cholesky_decompose,&
21 : cp_fm_cholesky_invert
22 : USE cp_fm_pool_types, ONLY: cp_fm_pool_p_type
23 : USE cp_fm_struct, ONLY: cp_fm_struct_create,&
24 : cp_fm_struct_release,&
25 : cp_fm_struct_type
26 : USE cp_fm_types, ONLY: cp_fm_create,&
27 : cp_fm_release,&
28 : cp_fm_type
29 : USE cp_log_handling, ONLY: cp_to_string
30 : USE cp_output_handling, ONLY: debug_print_level
31 : USE input_section_types, ONLY: section_vals_get,&
32 : section_vals_get_subs_vals,&
33 : section_vals_type,&
34 : section_vals_val_get,&
35 : section_vals_val_set
36 : USE kinds, ONLY: default_string_length
37 : USE message_passing, ONLY: mp_para_env_type
38 : USE optimize_basis_types, ONLY: basis_optimization_type,&
39 : flex_basis_type,&
40 : method_mo_fit_occ_virtual
41 : USE particle_types, ONLY: particle_type
42 : USE qs_density_matrices, ONLY: calculate_density_matrix
43 : USE qs_energy_init, ONLY: qs_energies_init
44 : USE qs_environment_types, ONLY: get_qs_env,&
45 : qs_environment_type
46 : USE qs_interactions, ONLY: init_interaction_radii
47 : USE qs_kind_types, ONLY: qs_kind_type
48 : USE qs_ks_methods, ONLY: qs_ks_update_qs_env
49 : USE qs_ks_types, ONLY: qs_ks_did_change
50 : USE qs_matrix_pools, ONLY: mpools_get
51 : USE qs_mo_io, ONLY: read_mo_set_from_restart
52 : USE qs_mo_types, ONLY: init_mo_set,&
53 : mo_set_type
54 : USE qs_rho_methods, ONLY: qs_rho_update_rho
55 : USE qs_rho_types, ONLY: qs_rho_get,&
56 : qs_rho_type
57 : USE string_utilities, ONLY: uppercase
58 : #include "./base/base_uses.f90"
59 :
60 : IMPLICIT NONE
61 : PRIVATE
62 :
63 : PUBLIC :: modify_input_settings, &
64 : allocate_mo_sets, &
65 : update_basis_set, &
66 : calculate_ks_matrix, &
67 : calculate_overlap_inverse
68 :
69 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'optbas_fenv_manipulation'
70 :
71 : CONTAINS
72 :
73 : ! **************************************************************************************************
74 : !> \brief change settings in the training input files to initialize
75 : !> all needed structures and adjust settings to basis optimization
76 : !> \param basis_optimization ...
77 : !> \param bas_id ...
78 : !> \param input_file ...
79 : !> \author Florian Schiffmann
80 : ! **************************************************************************************************
81 39 : SUBROUTINE modify_input_settings(basis_optimization, bas_id, input_file)
82 : TYPE(basis_optimization_type) :: basis_optimization
83 : INTEGER :: bas_id
84 : TYPE(section_vals_type), POINTER :: input_file
85 :
86 : CHARACTER(LEN=default_string_length) :: atom
87 : CHARACTER(LEN=default_string_length), &
88 13 : DIMENSION(:), POINTER :: abasinfo, obasinfo
89 : INTEGER :: ibasis, ikind, jkind, nbasis, nkind
90 : TYPE(section_vals_type), POINTER :: dft_section, feval_section, &
91 : kind_section, subsys_section
92 :
93 26 : feval_section => section_vals_get_subs_vals(input_file, "FORCE_EVAL")
94 13 : dft_section => section_vals_get_subs_vals(feval_section, "DFT")
95 13 : subsys_section => section_vals_get_subs_vals(feval_section, "SUBSYS")
96 13 : kind_section => section_vals_get_subs_vals(subsys_section, "KIND")
97 :
98 : CALL section_vals_val_set(feval_section, "PRINT%DISTRIBUTION%_SECTION_PARAMETERS_", &
99 13 : i_val=debug_print_level)
100 : CALL section_vals_val_set(dft_section, "SCF%PRINT%TOTAL_DENSITIES%_SECTION_PARAMETERS_", &
101 13 : i_val=debug_print_level)
102 : CALL section_vals_val_set(dft_section, "SCF%PRINT%DETAILED_ENERGY%_SECTION_PARAMETERS_", &
103 13 : i_val=debug_print_level)
104 :
105 : ! Add the basis containing the optimization templates to the basis-file list.
106 13 : CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", n_rep_val=nbasis)
107 13 : IF (basis_optimization%method == method_mo_fit_occ_virtual) THEN
108 : CALL section_vals_val_set(dft_section, "BASIS_SET_FILE_NAME", i_rep_val=nbasis + 1, &
109 4 : c_val=basis_optimization%template_basis_file)
110 : ELSE
111 : CALL section_vals_val_set(dft_section, "BASIS_SET_FILE_NAME", i_rep_val=nbasis + 1, &
112 9 : c_val=basis_optimization%work_basis_file)
113 : END IF
114 :
115 : ! Set the auxilarry basis in the kind sections
116 13 : CALL section_vals_get(kind_section, n_repetition=nkind)
117 36 : DO ikind = 1, nkind
118 : CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
119 23 : c_val=atom, i_rep_section=ikind)
120 23 : CALL uppercase(atom)
121 23 : CALL section_vals_val_get(kind_section, "BASIS_SET", n_rep_val=nbasis, i_rep_section=ikind)
122 23 : IF (nbasis > 1) THEN
123 : CALL cp_abort(__LOCATION__, &
124 0 : "Basis set optimization: Only one single BASIS_SET allowed per KIND in the reference input")
125 : END IF
126 : CALL section_vals_val_get(kind_section, "BASIS_SET", &
127 23 : c_vals=obasinfo, i_rep_val=1, i_rep_section=ikind)
128 23 : ALLOCATE (abasinfo(2))
129 23 : abasinfo(1) = "AUX_OPT"
130 23 : IF (SIZE(obasinfo) == 1) THEN
131 23 : abasinfo(2) = obasinfo(1)
132 : ELSE
133 0 : abasinfo(2) = obasinfo(2)
134 : END IF
135 : CALL section_vals_val_set(kind_section, "BASIS_SET", &
136 23 : c_vals_ptr=abasinfo, i_rep_val=2, i_rep_section=ikind)
137 23 : CALL section_vals_val_get(kind_section, "BASIS_SET", n_rep_val=nbasis, i_rep_section=ikind)
138 23 : CPASSERT(nbasis == 2)
139 :
140 92 : DO jkind = 1, basis_optimization%nkind
141 33 : IF (atom == basis_optimization%kind_basis(jkind)%element) THEN
142 :
143 23 : NULLIFY (abasinfo)
144 : CALL section_vals_val_get(kind_section, "BASIS_SET", &
145 23 : c_vals=abasinfo, i_rep_val=2, i_rep_section=ikind)
146 23 : ibasis = basis_optimization%combination(bas_id, jkind)
147 23 : CPASSERT(SIZE(abasinfo) == 2)
148 23 : CPASSERT(abasinfo(1) == "AUX_OPT")
149 23 : abasinfo(2) = TRIM(ADJUSTL(basis_optimization%kind_basis(jkind)%flex_basis(ibasis)%basis_name))
150 23 : EXIT
151 : END IF
152 : END DO
153 : END DO
154 :
155 13 : END SUBROUTINE modify_input_settings
156 :
157 : ! **************************************************************************************************
158 : !> \brief ...
159 : !> \param qs_env ...
160 : ! **************************************************************************************************
161 11 : SUBROUTINE allocate_mo_sets(qs_env)
162 : TYPE(qs_environment_type), POINTER :: qs_env
163 :
164 : INTEGER :: ispin
165 : INTEGER, DIMENSION(2) :: nelectron_spin
166 : LOGICAL :: natom_mismatch
167 11 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
168 11 : TYPE(cp_fm_pool_p_type), DIMENSION(:), POINTER :: ao_mo_fm_pools
169 : TYPE(dft_control_type), POINTER :: dft_control
170 11 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
171 : TYPE(mp_para_env_type), POINTER :: para_env
172 11 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
173 11 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
174 : TYPE(section_vals_type), POINTER :: dft_section
175 :
176 11 : NULLIFY (para_env)
177 : CALL get_qs_env(qs_env=qs_env, &
178 : dft_control=dft_control, &
179 : mos=mos, nelectron_spin=nelectron_spin, &
180 : atomic_kind_set=atomic_kind_set, &
181 : qs_kind_set=qs_kind_set, &
182 : particle_set=particle_set, &
183 11 : para_env=para_env)
184 11 : dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
185 :
186 11 : CALL mpools_get(qs_env%mpools, ao_mo_fm_pools=ao_mo_fm_pools)
187 22 : DO ispin = 1, dft_control%nspins
188 22 : IF (.NOT. ASSOCIATED(mos(ispin)%mo_coeff)) THEN
189 : CALL init_mo_set(mos(ispin), &
190 : fm_pool=ao_mo_fm_pools(ispin)%pool, &
191 11 : name="qs_env%mo"//TRIM(ADJUSTL(cp_to_string(ispin))))
192 : END IF
193 : END DO
194 :
195 : CALL read_mo_set_from_restart(mos, qs_kind_set, particle_set, para_env, id_nr=0, &
196 : multiplicity=dft_control%multiplicity, &
197 : dft_section=dft_section, &
198 11 : natom_mismatch=natom_mismatch)
199 :
200 11 : END SUBROUTINE allocate_mo_sets
201 :
202 : ! **************************************************************************************************
203 : !> \brief ...
204 : !> \param qs_env ...
205 : ! **************************************************************************************************
206 11 : SUBROUTINE calculate_ks_matrix(qs_env)
207 : TYPE(qs_environment_type), POINTER :: qs_env
208 :
209 : INTEGER :: ispin
210 11 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
211 : TYPE(dft_control_type), POINTER :: dft_control
212 : TYPE(qs_rho_type), POINTER :: rho
213 :
214 11 : NULLIFY (rho, dft_control, rho_ao)
215 :
216 11 : CALL qs_energies_init(qs_env, .FALSE.)
217 11 : CALL get_qs_env(qs_env, rho=rho, dft_control=dft_control)
218 11 : CALL qs_rho_get(rho, rho_ao=rho_ao)
219 22 : DO ispin = 1, dft_control%nspins
220 22 : CALL calculate_density_matrix(qs_env%mos(ispin), rho_ao(ispin)%matrix)
221 : END DO
222 11 : CALL qs_rho_update_rho(rho, qs_env)
223 11 : CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.TRUE.)
224 11 : qs_env%requires_mo_derivs = .FALSE.
225 11 : CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE.)
226 :
227 11 : END SUBROUTINE calculate_ks_matrix
228 :
229 : ! **************************************************************************************************
230 : !> \brief ...
231 : !> \param matrix_s ...
232 : !> \param matrix_s_inv ...
233 : !> \param para_env ...
234 : !> \param context ...
235 : ! **************************************************************************************************
236 27 : SUBROUTINE calculate_overlap_inverse(matrix_s, matrix_s_inv, para_env, context)
237 : TYPE(dbcsr_type), POINTER :: matrix_s
238 : TYPE(cp_fm_type), INTENT(OUT) :: matrix_s_inv
239 : TYPE(mp_para_env_type), POINTER :: para_env
240 : TYPE(cp_blacs_env_type), POINTER :: context
241 :
242 : INTEGER :: nao
243 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
244 : TYPE(cp_fm_type) :: work1
245 :
246 9 : CALL dbcsr_get_info(matrix_s, nfullrows_total=nao)
247 : CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=nao, ncol_global=nao, &
248 9 : para_env=para_env, context=context)
249 :
250 9 : CALL cp_fm_create(matrix_s_inv, matrix_struct=fm_struct_tmp)
251 9 : CALL cp_fm_create(work1, matrix_struct=fm_struct_tmp)
252 9 : CALL copy_dbcsr_to_fm(matrix_s, matrix_s_inv)
253 9 : CALL cp_fm_uplo_to_full(matrix_s_inv, work1)
254 9 : CALL cp_fm_cholesky_decompose(matrix_s_inv)
255 9 : CALL cp_fm_cholesky_invert(matrix_s_inv)
256 9 : CALL cp_fm_uplo_to_full(matrix_s_inv, work1)
257 9 : CALL cp_fm_struct_release(fm_struct_tmp)
258 9 : CALL cp_fm_release(work1)
259 :
260 9 : END SUBROUTINE calculate_overlap_inverse
261 :
262 : ! **************************************************************************************************
263 : !> \brief ...
264 : !> \param opt_bas ...
265 : !> \param bas_id ...
266 : !> \param basis_type ...
267 : !> \param qs_env ...
268 : ! **************************************************************************************************
269 644 : SUBROUTINE update_basis_set(opt_bas, bas_id, basis_type, qs_env)
270 : TYPE(basis_optimization_type) :: opt_bas
271 : INTEGER :: bas_id
272 : CHARACTER(*) :: basis_type
273 : TYPE(qs_environment_type), POINTER :: qs_env
274 :
275 : CHARACTER(default_string_length) :: elem
276 : INTEGER :: ibasis, ikind, jkind
277 644 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
278 : TYPE(dft_control_type), POINTER :: dft_control
279 : TYPE(gto_basis_set_type), POINTER :: gto_basis
280 644 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
281 :
282 : CALL get_qs_env(qs_env=qs_env, dft_control=dft_control, &
283 644 : atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set)
284 1869 : DO ikind = 1, SIZE(qs_kind_set)
285 4319 : DO jkind = 1, opt_bas%nkind
286 2450 : CALL get_atomic_kind(atomic_kind_set(ikind), name=elem)
287 2450 : CALL uppercase(elem)
288 3675 : IF (elem == opt_bas%kind_basis(jkind)%element) THEN
289 1225 : ibasis = opt_bas%combination(bas_id, jkind)
290 : CALL get_basis_from_container(qs_kind_set(ikind)%basis_sets, basis_set=gto_basis, &
291 1225 : basis_type=basis_type)
292 1225 : CALL transfer_data_to_gto(gto_basis, opt_bas%kind_basis(jkind)%flex_basis(ibasis))
293 1225 : CALL init_orb_basis_set(gto_basis)
294 : END IF
295 : END DO
296 : END DO
297 :
298 644 : CALL init_interaction_radii(dft_control%qs_control, qs_kind_set)
299 :
300 644 : END SUBROUTINE update_basis_set
301 :
302 : ! **************************************************************************************************
303 : !> \brief ...
304 : !> \param gto_basis ...
305 : !> \param basis ...
306 : ! **************************************************************************************************
307 1225 : SUBROUTINE transfer_data_to_gto(gto_basis, basis)
308 : TYPE(gto_basis_set_type), POINTER :: gto_basis
309 : TYPE(flex_basis_type) :: basis
310 :
311 : INTEGER :: ipgf, iset, ishell
312 :
313 2450 : DO iset = 1, basis%nsets
314 6337 : DO ishell = 1, basis%subset(iset)%ncon_tot
315 40071 : DO ipgf = 1, basis%subset(iset)%nexp
316 38846 : gto_basis%gcc(ipgf, ishell, iset) = basis%subset(iset)%coeff(ipgf, ishell)
317 : END DO
318 : END DO
319 10615 : DO ipgf = 1, basis%subset(iset)%nexp
320 9390 : gto_basis%zet(ipgf, iset) = basis%subset(iset)%exps(ipgf)
321 : END DO
322 : END DO
323 :
324 1225 : END SUBROUTINE transfer_data_to_gto
325 :
326 : END MODULE optbas_fenv_manipulation
|