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 Initialize a qs_env for kpoint calculations starting from a gamma point qs_env
10 : !> \par History
11 : !> 11.2016 created [JGH]
12 : !> \author JGH
13 : ! **************************************************************************************************
14 : MODULE qs_gamma2kp
15 : USE cp_control_types, ONLY: dft_control_type
16 : USE cp_dbcsr_api, ONLY: dbcsr_add,&
17 : dbcsr_p_type
18 : USE cp_subsys_types, ONLY: cp_subsys_type
19 : USE input_constants, ONLY: atomic_guess,&
20 : xc_none
21 : USE input_section_types, ONLY: section_vals_get_subs_vals,&
22 : section_vals_type,&
23 : section_vals_val_get,&
24 : section_vals_val_set
25 : USE kinds, ONLY: dp
26 : USE kpoint_types, ONLY: kpoint_create,&
27 : kpoint_type
28 : USE message_passing, ONLY: mp_para_env_type
29 : USE pw_methods, ONLY: pw_copy
30 : USE pw_types, ONLY: pw_c1d_gs_type,&
31 : pw_r3d_rs_type
32 : USE qs_energy_init, ONLY: qs_energies_init
33 : USE qs_environment, ONLY: qs_init
34 : USE qs_environment_types, ONLY: get_qs_env,&
35 : qs_env_create,&
36 : qs_environment_type,&
37 : set_qs_env
38 : USE qs_ks_methods, ONLY: qs_ks_update_qs_env
39 : USE qs_rho_atom_types, ONLY: rho_atom_type
40 : USE qs_rho_types, ONLY: qs_rho_get,&
41 : qs_rho_set,&
42 : qs_rho_type
43 : USE qs_scf_initialization, ONLY: qs_scf_env_init_basic
44 : USE qs_scf_types, ONLY: qs_scf_env_type
45 : USE scf_control_types, ONLY: scf_control_type
46 : #include "./base/base_uses.f90"
47 :
48 : IMPLICIT NONE
49 : PRIVATE
50 :
51 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_gamma2kp'
52 :
53 : PUBLIC :: create_kp_from_gamma
54 :
55 : ! **************************************************************************************************
56 :
57 : CONTAINS
58 :
59 : ! **************************************************************************************************
60 : !> \brief ...
61 : !> \param qs_env ...
62 : !> \param qs_env_kp ...
63 : !> \param with_xc_terms ...
64 : ! **************************************************************************************************
65 12 : SUBROUTINE create_kp_from_gamma(qs_env, qs_env_kp, with_xc_terms)
66 : TYPE(qs_environment_type), POINTER :: qs_env, qs_env_kp
67 : LOGICAL, OPTIONAL :: with_xc_terms
68 :
69 : INTEGER :: iatom, ispin, xc_func
70 : LOGICAL :: without_xc_terms
71 : TYPE(cp_subsys_type), POINTER :: cp_subsys
72 : TYPE(dft_control_type), POINTER :: dft_control
73 : TYPE(kpoint_type), POINTER :: kpoint
74 : TYPE(mp_para_env_type), POINTER :: para_env
75 : TYPE(qs_rho_type), POINTER :: rho, rho_gamma
76 : TYPE(qs_scf_env_type), POINTER :: scf_env
77 12 : TYPE(rho_atom_type), POINTER :: atoms(:), atoms_gamma(:)
78 : TYPE(scf_control_type), POINTER :: scf_control
79 : TYPE(section_vals_type), POINTER :: force_env_section, subsys_section, &
80 : xc_fun_section
81 :
82 : CALL get_qs_env(qs_env, &
83 : para_env=para_env, &
84 : input=force_env_section, &
85 12 : cp_subsys=cp_subsys)
86 :
87 12 : NULLIFY (subsys_section)
88 :
89 12 : NULLIFY (kpoint)
90 12 : CALL kpoint_create(kpoint)
91 12 : kpoint%kp_scheme = "GAMMA"
92 12 : kpoint%symmetry = .FALSE.
93 12 : kpoint%verbose = .FALSE.
94 12 : kpoint%full_grid = .TRUE.
95 12 : kpoint%eps_geo = 1.0e-6_dp
96 12 : kpoint%use_real_wfn = .TRUE.
97 12 : kpoint%parallel_group_size = 0
98 :
99 12 : without_xc_terms = .FALSE.
100 12 : IF (PRESENT(with_xc_terms)) without_xc_terms = .NOT. with_xc_terms
101 :
102 0 : IF (without_xc_terms) THEN
103 0 : xc_fun_section => section_vals_get_subs_vals(force_env_section, "DFT%XC%XC_FUNCTIONAL")
104 0 : CALL section_vals_val_get(xc_fun_section, "_SECTION_PARAMETERS_", i_val=xc_func)
105 0 : CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", i_val=xc_none)
106 : END IF
107 :
108 420 : ALLOCATE (qs_env_kp)
109 12 : CALL qs_env_create(qs_env_kp)
110 : CALL qs_init(qs_env_kp, para_env, cp_subsys=cp_subsys, kpoint_env=kpoint, &
111 : force_env_section=force_env_section, subsys_section=subsys_section, &
112 12 : use_motion_section=.FALSE., silent=.TRUE.)
113 :
114 12 : CALL get_qs_env(qs_env_kp, scf_control=scf_control)
115 12 : scf_control%density_guess = atomic_guess
116 12 : CALL qs_energies_init(qs_env_kp, calc_forces=.FALSE.)
117 :
118 12 : NULLIFY (scf_env)
119 12 : CALL qs_scf_env_init_basic(qs_env_kp, scf_env)
120 :
121 12 : CALL set_qs_env(qs_env_kp, scf_env=scf_env)
122 :
123 : ! Preserve the converged density fields, not the new environment's atomic guess.
124 12 : CALL get_qs_env(qs_env, rho=rho_gamma)
125 12 : CALL get_qs_env(qs_env_kp, rho=rho)
126 12 : CALL copy_frozen_density(rho_gamma, rho)
127 :
128 12 : CALL get_qs_env(qs_env, dft_control=dft_control)
129 12 : IF (dft_control%qs_control%gapw_xc) THEN
130 2 : CALL get_qs_env(qs_env, rho_xc=rho_gamma)
131 2 : CALL get_qs_env(qs_env_kp, rho_xc=rho)
132 2 : CALL copy_frozen_density(rho_gamma, rho)
133 : END IF
134 12 : IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) THEN
135 : ! The frozen GAPW potential also depends on these hard/soft onsite densities.
136 6 : CALL get_qs_env(qs_env, rho_atom_set=atoms_gamma)
137 6 : CALL get_qs_env(qs_env_kp, rho_atom_set=atoms)
138 6 : CPASSERT(SIZE(atoms) == SIZE(atoms_gamma))
139 20 : DO iatom = 1, SIZE(atoms_gamma)
140 14 : IF (.NOT. ASSOCIATED(atoms_gamma(iatom)%cpc_h)) CYCLE
141 14 : CPASSERT(ASSOCIATED(atoms(iatom)%cpc_h))
142 40 : DO ispin = 1, SIZE(atoms_gamma(iatom)%cpc_h)
143 14 : IF (.NOT. ASSOCIATED(atoms_gamma(iatom)%cpc_h(ispin)%r_coef)) CYCLE
144 14 : CPASSERT(ASSOCIATED(atoms(iatom)%cpc_h(ispin)%r_coef))
145 14 : CPASSERT(ASSOCIATED(atoms(iatom)%cpc_s(ispin)%r_coef))
146 11788 : atoms(iatom)%cpc_h(ispin)%r_coef = atoms_gamma(iatom)%cpc_h(ispin)%r_coef
147 11802 : atoms(iatom)%cpc_s(ispin)%r_coef = atoms_gamma(iatom)%cpc_s(ispin)%r_coef
148 : END DO
149 : END DO
150 : END IF
151 :
152 12 : CALL qs_ks_update_qs_env(qs_env_kp, print_active=.FALSE.)
153 :
154 12 : IF (without_xc_terms) THEN
155 : ! set back the functional
156 0 : CALL section_vals_val_set(xc_fun_section, "_SECTION_PARAMETERS_", i_val=xc_func)
157 : END IF
158 :
159 12 : END SUBROUTINE create_kp_from_gamma
160 :
161 : ! **************************************************************************************************
162 : !> \brief Copy converged Gamma density fields into the initialized k-point environment.
163 : !> \param source Converged Gamma density
164 : !> \param destination Initialized density with the same real/reciprocal grids
165 : !> \author Thomas D. Kuehne (tkuehne@cp2k.org)
166 : ! **************************************************************************************************
167 14 : SUBROUTINE copy_frozen_density(source, destination)
168 : TYPE(qs_rho_type), INTENT(IN) :: source
169 : TYPE(qs_rho_type), INTENT(INOUT) :: destination
170 :
171 : INTEGER :: idir, ispin
172 : LOGICAL :: drho_g_valid, drho_r_valid, rho_g_valid, &
173 : rho_r_valid, soft_valid, tau_g_valid, &
174 : tau_r_valid
175 14 : REAL(dp), DIMENSION(:), POINTER :: total_g_in, total_g_out, total_r_in, &
176 14 : total_r_out
177 14 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ao_in, ao_out
178 14 : TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_in, rho_g_out, tau_g_in, tau_g_out
179 14 : TYPE(pw_c1d_gs_type), DIMENSION(:, :), POINTER :: drho_g_in, drho_g_out
180 14 : TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_in, rho_r_out, tau_r_in, tau_r_out
181 14 : TYPE(pw_r3d_rs_type), DIMENSION(:, :), POINTER :: drho_r_in, drho_r_out
182 :
183 : CALL qs_rho_get(source, rho_ao_kp=ao_in, rho_r=rho_r_in, rho_g=rho_g_in, &
184 : tau_r=tau_r_in, tau_g=tau_g_in, drho_r=drho_r_in, drho_g=drho_g_in, &
185 : tot_rho_r=total_r_in, tot_rho_g=total_g_in, &
186 : rho_r_valid=rho_r_valid, rho_g_valid=rho_g_valid, &
187 : tau_r_valid=tau_r_valid, tau_g_valid=tau_g_valid, &
188 14 : drho_r_valid=drho_r_valid, drho_g_valid=drho_g_valid, soft_valid=soft_valid)
189 : CALL qs_rho_get(destination, rho_ao_kp=ao_out, rho_r=rho_r_out, rho_g=rho_g_out, &
190 : tau_r=tau_r_out, tau_g=tau_g_out, drho_r=drho_r_out, drho_g=drho_g_out, &
191 14 : tot_rho_r=total_r_out, tot_rho_g=total_g_out)
192 :
193 14 : CPASSERT(SIZE(rho_r_in) == SIZE(rho_r_out))
194 28 : DO ispin = 1, SIZE(rho_r_in)
195 14 : CALL pw_copy(rho_r_in(ispin), rho_r_out(ispin))
196 14 : CALL pw_copy(rho_g_in(ispin), rho_g_out(ispin))
197 : CALL dbcsr_add(matrix_a=ao_out(ispin, 1)%matrix, alpha_scalar=0.0_dp, &
198 14 : matrix_b=ao_in(ispin, 1)%matrix, beta_scalar=1.0_dp)
199 : ! XC fields need not be allocated when with_xc_terms is false.
200 14 : IF (ASSOCIATED(tau_r_out) .AND. tau_r_valid) THEN
201 6 : CPASSERT(ASSOCIATED(tau_r_in))
202 6 : CALL pw_copy(tau_r_in(ispin), tau_r_out(ispin))
203 : END IF
204 14 : IF (ASSOCIATED(tau_g_out) .AND. tau_g_valid) THEN
205 6 : CPASSERT(ASSOCIATED(tau_g_in))
206 6 : CALL pw_copy(tau_g_in(ispin), tau_g_out(ispin))
207 : END IF
208 14 : IF (ASSOCIATED(drho_r_out) .AND. drho_r_valid) THEN
209 0 : CPASSERT(ASSOCIATED(drho_r_in))
210 0 : DO idir = 1, SIZE(drho_r_out, 1)
211 0 : CALL pw_copy(drho_r_in(idir, ispin), drho_r_out(idir, ispin))
212 : END DO
213 : END IF
214 28 : IF (ASSOCIATED(drho_g_out) .AND. drho_g_valid) THEN
215 0 : CPASSERT(ASSOCIATED(drho_g_in))
216 0 : DO idir = 1, SIZE(drho_g_out, 1)
217 0 : CALL pw_copy(drho_g_in(idir, ispin), drho_g_out(idir, ispin))
218 : END DO
219 : END IF
220 : END DO
221 56 : IF (ASSOCIATED(total_r_in) .AND. ASSOCIATED(total_r_out)) total_r_out = total_r_in
222 14 : IF (ASSOCIATED(total_g_in) .AND. ASSOCIATED(total_g_out)) total_g_out = total_g_in
223 : CALL qs_rho_set(destination, rho_r_valid=rho_r_valid, rho_g_valid=rho_g_valid, &
224 : tau_r_valid=tau_r_valid .AND. ASSOCIATED(tau_r_out), &
225 : tau_g_valid=tau_g_valid .AND. ASSOCIATED(tau_g_out), &
226 : drho_r_valid=drho_r_valid .AND. ASSOCIATED(drho_r_out), &
227 58 : drho_g_valid=drho_g_valid .AND. ASSOCIATED(drho_g_out), soft_valid=soft_valid)
228 :
229 14 : END SUBROUTINE copy_frozen_density
230 :
231 : END MODULE qs_gamma2kp
|
